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ABSTltACT 

A  numerical  method  is  developed  based  on  potential  flow  non¬ 
linear  lifting  surface  theory  for  predicting  the  surface  velocities 
and  pressures  on  a  rotor  blade  of  an  arbitrary  helicopter  rotor  system 
which  is  executing  a  constant  rotational  and  constant  axial  translational 
motion  including,  specifically,  the  hover  flight  mode.  The  formula¬ 
tion  of  the  problem  is  exact  in  the  sense  that  the  normal  surface 
boundai^-  condition  is  satisfied  on  the  surface  of  the  rotor  blade. 

The  problem  is  governed  by  a  Fredholm  integral  equation  of  the  first 
kind  which  relates  a  singular  velocity  doublet  potential  surface 
distribution  applied  on  the  rotor  blades  and  wakes  to  the  normal  relative 
velocity  on  the  rotor  blade  surface.  The  wake  model  is  asaumea  to  be  of 
a  prescribed  shape. 

The  solution  of  the  integral  is  obtained  in  a  numerical  fashion 
by  approximating  the  actual  rotor  blade  upper  and  lower  surfaces  and 
wake  surfaces  by  a  finite  number  cf  elemental  surfaces  on  which  the 
doublet  strength  is  assumed  constemt  and  then  satisfying  the  resulting 
set  of  numerical  normal  surface  boundary  conditions  at  the  centroid  of 
each  of  the  blade  elemental  surfaces. 

A  computer  program  was  developed  for  the  lifting  surface  theory 
which  depends  on  a  given  geometry  of  the  wake.  I’liis  program  lends 
itself  to  an  iterative  procedure  for  a  future  force  free  wake  lifting 
surface  theox^^  analysis.  The  computer  results  for  two  case  studies  is 


also  presented. 
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which  they  sure  used.  The  physical  units  of  any  parameter  may  always 
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a)  All  lengths  are  in  units  of  feet. 

b)  All  forces  are  in  units  of  pounds. 
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by  the  rotor  rad: as  and  dividing  all  velocities  by  the  rotor  rotational 
tip  speed.  In  those  sections  where  it  becomes  necessary  to  distinguish 
dimensional  quantities  from  the  non-dimensional  quantities  we  have 
underlined  the  dimensional  quantity,  eg.  R  implies  units  of  I en^rth 
and  R  implies  a  non-dimensionalized  length. 

English 

A  Area 

B1  Hub  radius  teiken  to  be  the  distance  along  the  span 
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Chord  station  along  which  the  blade  spem  axis  lies. 

Chord  grid  station  expressed  in  terms  of  percent 
chord  length. 
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tip  speed  squared. 
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Angular  measure  that  two  adjacent  rotor  blade  span 
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Force 
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U  is  a  dummy  variable. 
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V  Vector  velocity  with  respect  to  the  (X,Y,Z,t)  reference 

frame. 
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to  the  (X,Y,Z,t)  reference  frame. 

W  Function  which  defines  wake  surface. 
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Indices 
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to  the  outboard  span  segment  first  on  the  upper  blade 
surface  then  on  the  lower  blade  surface 


I  Maalmum  number  of  indexed  blade  elemental  surfaces 

max 

per  blade. 

J  Same  as  I  above  and  used  in  terms  of  a  dxirainy 

summation  . 
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K  Indexed  wake  elemental  surface.  The  wake  elemental 
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Ktiay  Maximum  number  of  indexed  wake  elemental  surfaces 
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edge  span  segment . 

L  Indexed  blade  number.  The  blades  are  numbered  conse¬ 
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Number  of  rotor  blades  . 

M  Indexed  span  station  or  span  segment.  M=1  is  at  or 

near  the  root  section  respectively. 

Mfnny  Maximum  numoer  of  indexed  span  stations  . 

MI-l  Maximum  number  of  Indexed  span  segments.  (MM  =  M^^^  -  1) 
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surface  unit  normal  direction.  Also  used  to  index  the 
segment  of  the  elemental  surface  lying  between  corner 
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n  Maximum  number  of  i  indexed  corner  points  which  define 

an  elemental  surface  (usually  n  =  10 
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Reference  Coordinate  System 

(X,Y,Z)  Inertial  fixed  coordinate  system  in  general.  Used 

also  aa  a  dummy  reference  system  coincident  vfith  the 
(  /],  /J,  >'  )  systems  below 

(X' ,Y'  ,Z' )  Non-lnertial  body  fixed  coordinate  system 

(  /^  /  )  Element  fixed  coordinate  system 

(Xl,Yl,Zl)  Blade  one  body  fixed  coordinate  system  coincident  with 
(X' ,Y'  ,Z' ) 

i  )  Non-inertial  body  fixed  cylindrical  coordinate  system 
defined  in  the  usual  sense  within  (X',Y',Z') 
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CHAPTER  1.  INTRODUCTION 

1.1  Scope  and  Objectives 

There  is  a  need  for  a  more  rigorous  approach  to  the  treatment 
of  complex  three-dimensional  flows  for  geometries  of  certain  V/STOL 
aircraft.  In  the  case  of  rotors  In  quasi-hover,  the  downwash  velocity 
associated  with  the  generation  of  lift  is  large  when  compared  to  the 
axial  flight  velocity.  In  this  situation  the  classical  assimrptions 
of  lifting  potential  flow  aerodynamic  theoi*y  such  as  linearized  boundary 
conditions,  lifting  lines  and  rigid  non-force  free  wedces  do  not  lead 
to  accurate  predictions.  Without  accurate  inviscid  flow  predictions 
the  even  more  complicated  viscous  flow  analysis  cannot  even  be  begun. 

The  progress  in  high  speed  digital  computer  technology  now  allows 
one  to  formulate  the  flow  problem  more  realistically.  Although  the 
formulations  necessitate  the  approximation  of  the  integral  emd 
differential  equp.tlons  they  may  be  considered  exact  in  the  sense  that 
the  solution  is  atteiined  uniformly  as  the  computational  network  is 
refined. 

The  present  work  is  specifically  concerned  with  developing  a 
potential  flow  lifting  surface  theory  applicable  to  rotors  in  the 
axial  flight  mode  specifically  including  the  hover  mode.  The  theory 
is  an  exact  numerical  analysis  and  its  major  objective  is  to  predict 
the  local  three  dimensional  blade  surface  velocities  and  pressures. 

The  theory  necessarily  incorporates  a  prescribed  wake  model  because 
of  the  complex  nature  of  the  problem.  The  force  free  wake  analysis  is 
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to  be  achieved  by  successive  iterations  on  the  waXe  geometry  by 
incorporating  a  wake  prediction  method  into  the  analysis.  In  addition 
to  developing  the  actual  theory,  the  feasibility  of  applying  the 
theory  is  also  demonstrated  since  a  computer  program  was  also  developed 
and  is  presented  herein.  The  theory  and  program  developed  are  appli¬ 
cable  to  emy  eirbitrary  shaped  rotor  blade  having  a  finite  non- zero  hub 
radius  and  a  pointed  trailing  edge.  It  is  specifically  not  necessary 
for  the  rotor  blades  to  be  thin  as  the  surface  boundary  conditions  are 
satisfied  on  a  surface  network  described  on  the  wetted  blade  surface 
area.  Furthermore,  perturbation  velocities  etre  not  required  to  be 
small.  These  last  two  constraints  Eire  EiBsociated  with  the  so-called 
linearized  lifting  surface  theory  and  small  disturbance  theory 
respectively. 

1.2  Literature  Review 

In  the  discussion  of  any  lifting  surface  theory  one  must  first 
distinguish  between  two  basic  classes  of  problems.  One  clc'- j  of 
problems  is  concerned  with  the  prediction  of  local  surface  loadings 
by  assuming  a  loading  function  which  is  expressed  as  a  series  of 
assumed  modes  with  unknown  coefficients.  These  unknown  coefficients 
axe  then  obtained  by  satisfying  the  noimal  velocity  condition  either 
directly  or  indirectly  at  a  set  of  puints  whose  number  equals  the  number 
of  unknown  coefficients.  In  certain  cases  ^he  set  of  points  may  exceed 
the  number  of  unknown  coefficients  in  which  case  the  nonnal  velocity 
condition  is  satisfied  approximately  at  the  set  of  points  by 
appropriately  weighting  the  set  of  points.  Multhopp's  collocation 


method  (reference  1)  for  calculating  the  lift  distribution  of  wings 
in  subsonic  flows  exemplifies  this  method.  A  very  elegeint  analysis 
by  Verbaugh  (reference  2)  concerning  unsteady  lifting  surface  theory 
for  ship  screws  employs  the  acceleration  potentiad.  in  solving  this 
class  of  lifting  surface  theory  problem. 

The  second  class  of  lifting  surface  theory  problems  is  concerned 
with  the  prcf^lctlon  •'f  local  surface  loadings  by  assuming  a  distri¬ 
bution  of  surface  elements  on  which  the  loadings  are  unknown  but  oi. 
which  a  set  of  influence  coefficients  can  be  defined.  Perhaps  the 
most  complete  authoritative  discussion  of  this  method  is  that  presented 
by  Hess  and  Smith  (reference  3).  In  this  method  the  integral  equation 
resulting  from  the  application  of  the  normal  surface  boundeury  condition 
is  reduced  to  a  sum  of  integrations  to  be  performed  over  a  finite  set 
of  surface  elements  such  that  the  surface  boundary  condition  is  satis¬ 
fied  locally  at  one  point  on  each  surface  element.  The  loading 
function  is  some  velocity  potential  function  of  unknown  strength  which 
may,  however,  be  analytically  integrated  over  the  svirface  element 
region.  A  linear  set  of  equations  results  such  that  the  unknowns  are 
the  potential  function  strengths  on  each  surface  element  and  the 
coefficients  represent  the  elemental  integration  results.  The  results 
obtedned  by  this  method  are  excellent  as  documented  by  Hess  and  Smith 
preference  3)  for  nonlifting  bodies. 

The  lifting  surface  theory  methof^  presented  in  this  report  relies 
heavily  on  the  excellent  work  of  Hess  and  Smith  (reference  3).  The 
presented  problem  differs  fundamentally  in  that  this  theory  is  concerned 
with  a  rotating  lifting  body  behind  which  trails  a  wake  region. 


VaxlouB  authors  have  noted  the  complexity  of  the  calculations 
required  when  one  attempts  to  use  this  lifting  surface  theory  method 
fo.-  a  lifting  body.  As  late  as  1971  Johnson  (reference  1*)  noted  in 
substance  that  the  extent  of  the  calculations  involved  in  these  methods 
prohibited  the  direct  application  of  the  conventional  lifting  surface 
theory  technique  to  the  calculation  of  rotaiy  wing  air  loads.  Many 
authors  have  Ingeniously  attempted  to  tiimplify  the  exact  lifting  surface 
theory  method  in  order  to  attain  valid  results.  Erickson  (reference  5) 
reduced  (after  Prandtl)  the  lifting  blade  surface  model  to  a  lifting 
line  model.  As  such,  his  lifting  line  theory  was  baaed  on  a  bound 
vortex  line  and  a  continuous  wake  vortex  sheet  which  he  allowed  to 
distort  on  successive  iterative  steps.  The  contraction  pattern  was 
fixed  according  to  actuator  disc  theory.  Landgrebe  (reference  ^->)  alsc 
showed  that  the  readlstic  self-induced  distorted  wake  geometries  could 
be  computed  by  application  of  the  classical  Biot-Eavaxt  law  applied  over 
wake  vortex  filaments. 

Erikson  and  Hough  (reference  7)  showed  that  the  applicability  of 
the  lifting  line  model  for  hover  prediction  was  questionable  as  blade 
surface  induced  velocities  vary  rapidly  along  the  chord  direction 
which,  of  course,  would  invalidate  a  lifting  line  model.  The  reason  for 
this  rapid  variation  lies  in  the  fact  that  the  wake  has  a  pronounced 
influence  on  the  rotor  blade  because  of  its  near  proximity  in  hover. 

At  Sikorsky  Aircraft  Rorke  and  Wells  (r<»^erence  8)  have  described 


another  unique  variation  on  the  true  rotor  lifting  surface  theory. 

They  have  coupled  a  prescribed  wake-momentum  analysis  to  the  conventional 
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Btrip-momentum  thaory  in  order  to  predict  the  rotor  hover  performance. 

The  prescribed  waJce  geometry  in  this  method  ia  determined  in  part 
by  a  theoretical  analysif ,  the  details  of  which  were  presc^nted  by 
Clark  and  helper  (rejr>'jnce  9)*  This  analysis  is  a  true  engineering 
design  euialysis  emd  has  been  optimized  so  as  to  require  very  little 
computer  time.  This  technique  does,  of  course,  require  airfoil 
sectional  aerodynamic  coefficients. 

There  are  other  variations  )f  the  rotor  lifting  surface  theory 
presented  in  the  literature  but  to  this  authors  knowledge  none  of  the 
so  caUed  rotor  lifting  surface  theories  presented  are  in  fact  true 
applications  of  the  ideal  lifting  surface  model.  Furthermore,  it 
appears  that  no  single  reported  rotor  prediction  method  is  capable  of 
piedicting  local  surface  velocities  or  pressures  on  some  erbitrary 
rotor  geometry  surface.  Thus  it  appears  that  design  studios  of  new 
rotor  blade  shapes  differing  significantly  from  existing  blade  shapes 
cannot  be  performed  at  the  present  time  with  any  level  of  confidence. 
Because  of  this  technological  dej.ioit  (see  references  10  and  11) 
it  was  decided  to  attempt  to  develop  a  true  rotor  lifting  surface 
theory  and  actually  apply  this  theory  in  terras  of  an  exact  numerical 
sense.  This  work  is  concerned  with  the  initial  phase  of  the  development, 
that  is,  for  a  prescribed  wed^e  trailing  an  arbitrary  shaped  body 
develop  a  lifting  surface  theory  which  will  predict  for  the  axial  flight 
mode  local  surface  velocities  and  loadings.  The  succeeding  phase  will 
be  to  use  the  theory  and  program  of  the  initlil  phase  and  modify  them 
so  as  to  include  a  wake  iterative  scheme  in  order  to  include  a  force 
free  wedte  aneilysis  into  the  lifting  surface  theory. 


L.-,. 


We  have  not  attempted  to  review  here  the  subject  of  lifting 
surface  theory  in  its  entirety  but  have  rather  restricted  the  review 
to  selected  current  rotor  lifting  surface  theories  indicative  of  the 
general  development  trends.  For  an  additional  literature  review 
concerned  mainly  with  lifting  surface  theory  applied  to  planar  flows 
refer  to  DJoJodihardJo  (reference  12).  DJoJodihardJo  and  Widnall 
(reference  13)  in  a  paper  which  sximmarizes  the  previous  reference 
presents,  in  part,  a  discussion  of  the  doublet  velocity  potential 
which  we  have  used  to  verify  the  derived  integrated  doublet  velocity 
potential  used  herein. 

In  addition  to  the  explicit  references  above  w  have  included  a 
list  of  references  which  we  have  used  for  obtaining  fundamental 
information  and  for  obtaining  information  related  to  general  rotor 
performance  prediction  methods. 

1 . 3  Report  Layout 

We  have  presented  in  Chapter  2  the  formulation  of  the  problem  in 
terms  of  the  governing  equations.  Chapter  3  presents  a  step  by  step 
discussion  of  the  overall  problem  solution.  Chapter  discusses  the 
computer  program  in  a  general  manner.  The  results  of  two  computer  run 
caaes  are  presented  in  Chapter  5.  In  Cnapter  6  we  have  discussed 
extensions  to  the  present  problem.  We  have  relegated  all  discussion 
material  not  actually  essential  to  the  main  problem  discussion  to 
Appendices  A,  B  and  C  so  as  not  to  interrupt  the  overall  problem 
discussion.  Appendix  D  describes  the  computer  program  in  detail  as  to  its 
options  and  input/output. 
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CHAPTER  2.  DISCUSSION  OF  THE  BASIC  PROBLEM 

2.1  Formulation  of  the  Governing  Equations 

As  a  lifting  rotor  moves  into  the  air  it  disturbs  the  air  in 
such  a  manner  as  to  derive  its  lift.  This  problem  is  concerned 
with  the  prediction  of  the  local  surface  pressure  acting  on  the 
rotor  blades  in  hover  or  axial  flight  ti  rough  an  analysis  of  the 
rotor  induced  velocities.  We  shad.!  formuicte  in  this  section  the 
equations  governing  the  fluid  motion. 

Let  us  consider  the  lifting  rotor  system  to  consist  of: 

a)  a  three  dlmonsional  body  of  arbitrary  shape  which  is 
executing  a  constant  rotatory  eind  translatory  motion,  and 

b)  a  wake  which  trails  the  lifting  body.  The  surface  of  the 
body  may  be  represented  by 

1^)  ~  O  (2.1.1) 

and  ihe  wake  following  the  body  may  be  defined  by  a  surface  of  velocity 
discontinuity  given  by 

lA/  ^  tj  ~  O  (2.1.2) 

In  the  above  equations  R  is  the  position  vector  with  respect  to 
some  fixed  Inertial  reference  frame.  The  external  flow  field  is  assumed 
to  be  an  incompressible  flow  field  which  is  inviacid  and  initially 
irrotatlonal  and  at  rest.  The  wake  is  further  assumed  to  be  composed 
of  two  surfaces  coincident  with  each  other.  Each  weike  surface  is 
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assigned  to  have  its  origin  at  some  Infinitesimal  region  located  on 
the  body  upper  and  lover  surfaces  at  the  trailing  edge.  These  trailing 
edge  regions  as  such  represent  the  lines  along  which  the  viscid 
boundary  layer  smoothly  leaves  the  trailing  edge.  It  is  assumed  in 
the  analysis  to  follow  that  there  exists  no  flow  separation  from  the 
body  except  at  the  trailing  edge,  thus  the  body  must  have  a  sharp 
trailing  edge. 

It  follows  row  from  the  condition  of  irrotationality  .  d  the 
continuity  equation  that  a  velocity  potential  can 

be  defined  such  that  it  must  satisfy  Laplace's  equation 

^  .  (2.1.3) 

The  velocity  potential  at  aui  exterior  field  point  can  be  given  by  an 
Integral  equation  which  incorporates  a  distribution  of  singularities  over 
the  lifting  body  and  wake  surface.  In  this  analysis  we  will  take  as 
our  distributed  siurface  singularities  the  doublet  or  dipole  which  is 
itself  composed  of  two  more  basic  singul.-rlties ,  namely  the  source 
and  sink.  In  Appendix  C  is  presented  a  discussion  of  the  doublet 
potential  and  its  axis  convention.  The  resulting  integral  equation 
for  the  doublet  velocity  potential  is  given  by 


fMyJ-- 


(2.1. I) 
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where 

/?  is  the  position  vector  to  some  field  point,  P. 
is  the  position  vector  to  some  surface  point. 

is  the  unit  outward  surface  normal. 

■  X’  -  Rs  »  which  is  the  vector  from  the  surface 
point  to  the  field  point. 

is  the  doublet  strength  per  unit  area  at  some  surface 
point . 

The  potential  as  given  satisfies  Laplace's  equation  identically. 
The  above  velocity  potential  is  subject  to  the  following  boundary 
conditions : 

a)  In  the  far  region  away  from  the  doublet  surface  distribution 
the  fluid  velocity  (V)  given  by 

^  (2.1.5) 

should  approach  zero.  Thus  the  fax  boundary  condition  becomes 

^  O  .  (2.1.6) 

This  boundary  condition  is  inherently  satisfied  by  the  doublet 
velocity  potential. 

b)  In  the  region  of  the  body  surface  the  normal  velocity  at  the 
surface  must  be  zero.  This  kinematic  normal  boundary  condition 
may  be  expressed  as 
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dt 


-  O  . 

(2.1.7) 


c)  In  the  region  of  the  valce  surface  -  since  the  wake  cannot 
maintain  a  pressure  discontinuity,  the  pressure  across 
the  wake  surface  must  be  continuous .  Thus 


^  ^  ^  (2.1.8) 

on 

where  and  Pj,  sure  the  pressures  on  the  local  wake 

upper  and  lower  surfaces. 

d)  In  the  region  of  the  trailing  edge  the  fluid  velocity  must 
be  directed  smoothly  from  the  body  surface  to  the  wake 
surface.  This  is  the  Kutta-Joukowskl  condition. 


The  pressure 
equation  of  motion 


P  anywhere  in  the  flow  field  is  given  by  the 


VP 


(2.1.9) 


which  is  uncoupled  from  the  governing  kinematic  equation  given  by 
Laplace's  equation.  It  is,  however,  coupled  into  the  overall  problem 
through  the  boundary  condition  (c)  above. 

For  convenience  we  will  now  transform  the  above  equations  which 
are  expressed  with  respect  to  a  fixed  inertial  reference  freune  to  a 
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non-lnertlal  body  reference  frame  designated  by  primed  coordinates. 
Let  us  assume  that  the  body  Is  treinslatlng  and  rotating  about  the 
translational  axis  with  respect  to  the  Inertial  reference  frame  at 
a  constant  translational  velocity  given  by  Vj-  and  consteint 
rotational  velocity  given  by  ^  .  In  Appendix  A  Is  presented  a 

summary  of  the  trauisformation  relations.  The  body  surface  and  wake 
surface  equations  given  by  equations  (2.1.1)  and  (2.1.2)  transiorm 
respectively  to 


(2.1.10) 

(2.1.11) 


The  velocity  potential  In  the  body  fixed  reference  frame  becomes 
independent  of  time  and  Is  given  by 


Laplace's  equation  transforms  directly  to 


(2.1.13) 


Whereas  before,  however,  was  the  total  fluid  velocity 

V  »  »  in  the  body  fixed  system  becomes  the 

velocity  induced  about  a  relative  free  stream  velocity.  If  we  define 


v'^  //^'J 


(2.1.1U) 


aB  the  relative  free  stream  velocity  given  by 


and  1/ 


-  I  1<J  jK  M  '  ^  J 


(2.1.15) 


then  the  Inertial  and  body  fixed  velocities  are  related  as  follows 


)/  -  i/  ^  ^  M  /■  U, 


(.^.1.16) 


where 


_  / 
i/ 


(2.1.17) 


The  integral  equation  for  the  doublet  velocity  potential  becomes 


(2.1.18) 


where  the  symbols  are  the  same  as  before  except  the  reference  coordinate 
system  is  now  the  primed  body  fixed  reference  system.  The  boundary 
conditions  in  the  new  reference  plane  become: 
a)  In  the  far  region 


(2.1.19) 


which  remains  Inherently  satisfied. 


b)  On  the  body  atirface  since 


L  Z  i. 


(2.1.20) 


and  since 
relation 


JT 


ve  must  satisfy  the  following 


-  -  1/  '  V5 

S/  A'  J  ~  O 

(2.1.21) 


c)  On  the  wake  surface. 


(2.1.22) 


The  Kutta  condition  (boundary  condition  d)  may  be  applied  as  stated 
earlier.  The  pressure  anywhere  in  the  flow  field  may  be  determined  from 
the  transformed  equation  of  motion  which  will  be  discussed  later. 

Since  the  velocity  potential  is  a  solution  to  Laplace's  equation 
and  since  it  also  inherently  satisfies  the  far  boundary  condition  we 
no  longer  need  consider  these  two  equations.  Let  us  consider  now  the 
near  boimdary  condition.  If  we  define  to  be  the  outward  surface 

unit  normal  vector,  then  we  may  write 


Purthennore,  using  this  result  and  equation  (2.1.1U)  ve  may  revrite 
the  near  boundary  condition  given  by  equation (2. 1.21) as  follows: 

y  ^  ^  (2.1.2i4) 

on  -  6>  . 

Written  in  this  form  we  see  now  that  the  solution  to  the  given  proble’i 
requires  solving  in  the  body  fixed  reference  system  an  integral 
equation  given  by  the  near  boundary  condition  expressed  by  equation 
(2.1.2^)  where  is  defined  by  equation  (2.1.18)  subject  to  explicit 

boundary  conditions  (c)  and  (d).  Note  that  the  solution  does  not 
necessarily  involve  the  potential  but  rather  involves  the  vector 

gradient  of  ^  ,  that  is  V  .  In  essence  then  the  formulation 

of  the  basic  governing  equations  is  complete. 

Before  proceeding  into  the  detailed  formulation  we  shall  first 
discuss  the  doublet  potential  in  order  to  further  clarify  its  use. 
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2.2  Doublet  Potential  Dlacueslon 


We  have  taken  as  our  potential  function  a  potential  derived  on 
the  basis  of  a  surface  doublet  (dipole)  distribution.  In  particuleur  ve 
note  using  Figure  1  that  the  potential  Induced  at  a  point  P  by  a  surface 
doublet  distribution  given  by  whose  axis  is  everywhere  normal 

to  the  local  surface,  that  is 


(2.2.1) 


is  given  by 


(2.2.2) 


In  Appendix  C  we  have  derived  this  doublet  potential  as  presented 
in  terms  of  its  more  basic  source  plus  sink  potentials  in  order  to 
clarify  the  axis  convention  of  the  doublet. 

We  have  chosen  to  model  the  rotor  and  wake  surfaces  with  a  doublet 
distribwiLion  rather  than  a  vorticlty  distribution,  which  could  also 
have  been  used,  for  a  number  of  reasons.  The  main  reason  stems  from  the 
fact  that  a  vorticity  model  does  require  the  specification  of  two 
functions  lying  along  the  surface  in  two  vector  directions  mutually 
perpendicular  to  each  other.  This  aspect  of  vorticity  in  light  of  the 
doublet  model  would  unnecessarily  complicate  the  geometry  of  the  problem. 
The  two  models  are  of  course  mathematically  related  (refer  to  Appendix 
B)  and  both  inherently  satisfy  Laplace's  equation  and  the  boundary 
condition  at  infinity. 


lf> 


The  application  of  the  doublet  model  Involves  a  doublet  distribution 
applied  on  the  wetted  areas  of  all  rotor  blades.  Furthermore,  trailing 
continuously  from  the  upper  €md  lower  blade  surfaces  are  two  sheets 
'  pon  which  are  also  applied  doublets.  These  two  sheets  are  taken  to 
be  coincident  with  each  other  thus  allowing  for  vector  suras  of  the 
doublet  strengths. 

It  should  be  noted  at  this  point  that  the  doublets  are  distributed 
on  the  blade  surfaces.  This  method  is  to  be  distinguished  from  the 
so  ceilled  linearized  method  which  through  linearization  of  the  surface 
boundeury  condition  would  allow  doublets  to  be  distributed  on  a  mean 
camber  plajie.  The  blade  surface  distribution  although  more  complex  is 
of  course  exact  and  does  not  Impose  the  planform  restrictions  of  the 
linearized  method. 
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2.3  Evaluation  of  the  Doublet  Wake  Strength 

Let  us  consider  now  the  doublet  velocity  potential  in  more  detail. 
Suppose  ve  first  separate  the  wake  region  integration  from  the  blade 
region  integi.‘ation  such  that  ve  may  write 


(2.3.1) 


where  now  represent  the  doublet  strengths 

on  the  body  and  wake  regions  respectively.  As  was  pointed  out  earlier 
the  wake  is  actually  composed  of  two  coincident  sheets  on  each  of  which 
the  doublet  strength  is  given  as  an  upper  doublet  strength 
auid  a  lover  doublet  strength  .  If  we  account  for  the 

opposite  direction  of  the  doublet  cucis  and  take  as  a  net  eucis  direction 
the  outward  normal  of  the  upper  wake  surface  then  in  effect  the  net 
doublet  strength  on  the  wake  surface  is  given  locally  by 


•  (2.3.2) 

In  Appendix  B  it  is  shown  that  Helmholtz’s  theorem  on  conservation  of 
vortlclty  logically  extends  to  the  conservation  of  wake  doublet  strength. 
In  this  problem  the  blade  loadings  are  Independent  of  the  rotational 
azimuth  body  position  and  as  such  the  wake  doublet  strength  remains 
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constant  edong  a  wake  streamline.  Following  along  a  wake  streamline 
up  to  the  trailing  edge  by  utilizing  the  Kutta  condition  we  find  that 
along  any  wake  streamline  that  originates  at  trailing  edge  span 
station  the  doublet  wake  strength  is  constant  and  is  given 

by  the  difference  in  value  of  the  blade  doublet  upper  and  lower  trailing 
edge  surface  doublet  strength,  that  is 

=  -^6  .  (J.3.3) 

Thus  the  wake  doublet  strength  varies  with  spem  station,  is  constant 
along  a  streamline  and  is  given  in  terms  of  the  body  trailing  edge 
doublet  strengths.  The  velocity  potential  may  now  be  written  as 


(2.3.»4) 
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Nuiurlcal  Reduction  of  the  Integral  Expreaelon  for  the  Doublet 
PotentliJ.  and  Ita  Darivatlve 

We  may  aov  reduce  this  integral  expression  for  the  doublet 
potential  to  a  numerical  relation.  In  order  to  do  this  we  first  eissume 
that  these  ere  equally  spaoed,  identical  rotor  blade  bodies  whose 
surface  equations  are  given  by 

/a'J  ==  O  Z-'J,  .  (2.U.1) 

SlmlleLTly  we  assume  that  there  are  Lmny  equally  spaced,  identical  wake 
siirfaces  whose  surface  equations  cure  given  by 

^  Z?  .  (2.1*. 2) 

T'hus  equation  (2.3.1^)  may  now  be  written  as 


(2.1*. 3) 
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Let  us  assume  now  that  rotor  blade  body  one  and  wake  surface  one  can 
be  subdivided  into  M  =»  1,  MM  span  segments.  Furthermore  we  shall 
assume  that  along  blade  segment  M  there  are  N  «  1,  NNjj  and  NNl 
upper  and  lower  surface  planar  areas  respectively.  In  the  wake  one 
region  we  assume  that  along  wake  span  station  M  there  are  K  =  1 , 

Ktnny  streamline  defined  planar  wal:e  surface  areas.  We  assume  further 
that  on  blade  bodies  L  =  2,  Ljnax  wake  surfaces  L  =  2,  Lmax  there 
are  similarly  defined  elements  which  are  determined  from  blade  body 
one  and  wake  surface  one  by  simple  element  rotation.  If  in  addition 
we  assume  that  the  doublet  strength  on  any  element  is  a  constant  we 
may  rewrite  the  doublet  potential  given  by  equation  (2. >4.3)  as  follows 


-- 


- '  I  -  / 


M 


5/^ 


///V 


yr  ^  1  —  f 


noting  now  that  the  blade  loadings  and  thus  the  doublet  strengths  are 
independent  of  blade  azimuth  position  in  this  problem  we  may  identify 


^4 


--  ^^6 


= 


Lt  1 


Z/',  /^J. 


(2.14.5) 


igiunurniamiKm 


PI 


Thp.t  is  the  doublet  strength  on  blade  L  upper  surface  element  (N,M) 
has  the  seme  value  of  doublet  strenp.th  on  blade  L+1  upper  surface 
element  (N,M).  Applying  this  condition  also  to  the  lower  blade  surface 
elements  and  to  the  wake  elements  we  find  that  the  velocity  potential 
may  now  be  written  as 
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The  induced  velocity  at  some  field  point  alonp  a  direction 
may  be  written  as 


(2.I4.O) 


where  the  derivative  with  respect  to  S  of  the  potential  follows  directly 


from  equation  (2.14.6)  and  is  given  by 


A//'/ 

■  A  -r  A  ^ 


^  ^y'y 


/  -  /  X  - 


(2.I4.9) 


where  we  have  defined 


(2.U.10) 


Here  J^yOt  ,  is  a  velocity  influence  coefficient  in  the  S 
'y- 

direction  whi^  represents  the  geometric  influence  of  some  surface 
element(l,J)  of  blade  or  wake  L  on  some  point 


’^■ 
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2.5  Numerical  Reductidn  of  the  Integral  Surface  Boundary  Condition 

Let  US  summarize  at  this  point  the  major  results  of  the  previous 
discussions.  First  it  has  been  shown  that  the  equation  to  be  solved 
is  the  near  boundary  condition  given  by  equation  (2. 1.2*4)  and  repeated 
here  is 

^  (2.5.1) 

^/y  S^/'J  ^ 

This  equation  states  that  the  normal  fluid  velocity  on  the  body  must 
be  zero.  In  addition  this  equation  is  subject  to  the  following 
additional  constraints: 

A)  The  fluid  flow  must  leave  the  trailing  edge  smoothly. 

B)  The  pressure  must  be  continuous  across  the  wake.  This 
ensures  the  exiotance  of  a,  so  called,  force  free  wake. 

In  section  2.h  it  was  shown  that  the  integral  expression  for 
and  could  be  reduced  to  a  numerical  equation  involving  the 

summation  of  finite  elemental  surface  contributions.  On  each  of  these 
elemental  surfaces  it  was  assumed  that  the  local  doublet  strength  was 
a  constant. 

Let  us  now  relax  the  boundary  condition  expressed  by  equation 
(2.5.1)  such  that  instead  of  requiring  it  to  hold  on  the  entire  body 
region  ^ we  now  require  it  to  hold  only  at  a  finite 
number  of  body  points.  We  will  further  take  these  control  points  to 
be  located  at  say  the  centroid  of  the  previously  defined  elemental 
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body  areas .  Because  of  the  symmetry  involved  we  need  only  require 
the  normal  velocity  to  be  zero  on  one  blade's  control  points  in  order 
to  ensure  a  zero  velocity  on  all  other  corraspondint;  blade  control  points.  In 
previous  discussions  for  preciseness  we  identified  the  upper  and  lower 
surface  elemental  areas  on  a  perticular  blade  as  element  ^ 

Let  us  now  identify  these  same  blade  elements  here  using  a  single 
index  say  1  =  1,  .  Thus  the  numerical  form  of  the  near  boundary 

condition  becomes  a  set  of  ImAy  equations  such  that 


(2.5.2) 

for  all  I  =  1,  Iniax 


where 


represents  the  vector  gradient  of  the  potential 
on  the  centroid  of  element  I. 


represents  the  free  streeun  velocity  on  the  centroid 


of  element  I 


represents  the  outward  unit  normal  of  element  I 
at  its  centroid. 

_  /  — 

Since  ^ ^  Just  the  induced  velocity  along  the  normal 

direction  of  element  I  evaluated  at  the  centroid  of  element  I  we  may 
write  that 


(2.5.3) 


Now  using  equation  (2.U.9)  we  may  evaluate 


LfP,;P^J.  i.l*  IWJ,«lAi.l,IHill||4l!iUill|, 
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such  that 


we  have 


-ji  2  /■ 

V/^  Jr]  V 


/y/y 
/y--/ 


'•i/ 


ijr 


/r--/ 


/V. 


(2.5.i^) 


where  the  first  term  represents  the  body  surface  induced  normal  velocity 
on  control  point  I  and  the  second  term  represents  the  wake  induced 
normal  velocity  on  control  point  I.  Note  that  the  normal  influence 
coefficients  are  a  function  of  the  geometry  of  the  problem  only.  Note 
further  that  singularities  in  the  evaluation  of  particular  influence 
coefficients  will  occur  when  we  attempt  to  evaluate  an  elemental 
influence  on  itself.  The  actual  evaluation  of  these  influence 
coefficients  and  the  resolution  of  the  singularities  will  be  presented 
in  section  2.6. 

If  we  now  substitute  this  equation  into  equation  (2.5.2)  using 
equation  (2.5.3)  we  may  write  the  final  form  of  the  1  =  1,  Imax  set  of 
near  surface  boundary  conditions  to  be  satisfied  as  a  numerical  set  of 


equations  given  by 


/V/7 

'  ^  /V.  -  ^  7^  ^ 

/y-/  ^*^  A--J  /‘^ 


/  ^  ^yy,  ,  -  V 


"  A 

J 


/Aw  >^  '  -^  , 


(2.5.5) 


Given  now  the  free  stream  conditions  in  terms  of  the  rotor  climbing 

speed  and  rotor  rotational  rate,  and  given  the  blade  geometry  and  some 

prescribed  wake  geometry  then  the  above  set  of  equations  constitute  a 

set  of  Imax  lineeu:  independent  eilgebraic  equations  for  the  Imax 

unknown  blade  element  doublet  strengths  given  by  above.  Uote 

that  the  wake  doublet  strengths  are  given  in  term?  of  the  blade  upper 

and  lower  surface  trailing  edge  doublet  strengths  and  as  such  do  not 

constitute  any  additional  unknowns  as  they  may  be  expressed  in  terms  of 

the  blade  element  doublet  strengths  ^  .  Once  the  solution 

values  in  terms  of  are  obtained  then  the  velocities  at  any 

'J 

point  'J  in  any  specific  direction  may  be  obtained  using 
equation  {2.k.9).  The  potential  likewise  may  be  obtained  using  equation 
(2.1*.6).  Given  the  velocities  at  a  point  then  the  pressure  may  be 
calculated  as  is  shown  in  section  2.7-  In  essence  then,  this  is  the 
problem  which  we  have  solved.  Note  that  we  ha’  .*  incorporated  the 
tangency  flow  constraint  into  this  problem  directly.  Neither  have 
we  enforced  the  force  free  wake  constraint  into  the  prcoiem.  These 
constraints  are  to  be  satisfied  (in  subsequent  work)  by  an  iterative 
technique  which  iterates  on  the  wake  geometry  based  on  a  previous 
solution  whose  wake  geometry  was  prescribed.  If  the  prescribed  wake 
is  coincident  with  the  real  wake  then  the  tangency  flow  constraint  and 
ajid  the  zero  pressure  discontinuity  across  the  wake  SLre  inherently  satisfied. 
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In  the  preceding  discussion  we  selected  the  control  points  at 
which  the  normal  surface  boundary  conditions  are  to  be  satisfied  as 
the  centroids  of  the  elemental  surfaces.  We  actually  chose  the  cen¬ 
troid  location  for  two  reasons.  Firstly,  the  centroid  of  an  elemental 
planar  area  is  most  logically  the  point  which  best  represents  the 
surface  area  and  secondly,  this  point  is  the  point  on  an  elemental 
planar  surface  where  the  self  induced  doublet  normal  velocity  magnitude 
is  a  minimum  and,  as  such,  the  various  elemental  centroids  represent 
consistent,  comparable  points. 
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follows 
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^  (2.6.3) 

where 

JZ/’J  -  (2.6. 14) 


-  y^-'/-‘  /- /^ 'I i  '^'-J 

(2.6.5) 

\J=  /y  j  j  /y  j-A/j  7/>'y^  v/  / 


(2.6.6) 

- 

-1: 

(2.6.7) 

(2.6.8) 

- 

ZP  -  P,  J/'/  y*/,  J 

(2.6.9) 

/ jx- /.  y  J'y  y/.  J' 

(2.6.10) 

A'A*  A  yry/vy  yy^j/<£>//yy 

y / y> i' A  y 

A'yy/f^^^yy  >  ty/ynyA/j  zp/./'/^^, 
y//y  A/y4y?r/7y  py/v/r  y/yzAf^/L . 
yy  /.  --  />  7/yy  yz  =  Z.  (2.6.11) 

(2.6.12) 

If  the  field  point  P  lies  on  the  plane  of  the  surface  doublet  distri¬ 


bution.  that  is  when  ^  -  9/.  -  ^  ,  then  the  value  of  the  integral 

taken  as  the  limiting  process  as  Jr  approaches  Z  positively 


results  in 


tfiiiiiliiiiilliiljiii 


Httittiuiiiliiii 


. . 


^/^J  -  -  /  For  F  within  doublet  surface  (2.6.13) 

^//>J  -  /7  Por  P  oatniJe  doublet  surface  (2.6. ll^) 

For  P  on  doublet  surface  edge.  (2.6.15) 

Given  this  expression  for  we  can  now  evaluate 

as 

V-fZ/’J  -  -  .  (j.6.16) 


This  operation  is  ceirried  out  in  Appendix  C.  The  results  when 
referred  to  cartesian  coordinates  are 

~  '  ~  ^7/  (2.6. IT) 

where  ^  is  a  dummy  independent  variable  standinp;  for  x,  y  or  7,  and 

/f 

yJ/zdV  />v  (2.6.18) 

7‘  -A/jJ-  /x^/rv  A  '  Aj  J i  Z/Vy^/^  y 

(2.6.19) 

(2.6.20) 


A  i  Aj)^  /  A x^  j-  yj^  ~  Ay  J  ^ 


--  ^  --  h  *''  ■  A,  }>■  7 


where 

Aaaj  = 
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The  symbols  used  here  are  the  seune  as  given  following  equation  (2.6.3). 
The  value  for  is  to  be  taken  as  zero  in  the  event 

The  relation  for  Jfj  /f'Y is  valid  for  all  field  points  whether 
they  lie  on  or  off  the  plane  of  the  doublet  distribution  except  for 
the  case  where  P  lies  on  the  edge,  in  which  case,  the  integred  is 
undefined.  Note  also  that  regardless  of  the  size  or  shape  of  the  planar 
doublet  distribution  equation  (2.6.17)  shows  that 


(2.6.21) 


M  ^  (2.6.22) 


The  normal  derivative  with  respect  to  the  doublet  surface  is  well 
behaved  everywhere  except  for  the  edge  where  it  is  undefined.  It  is 
continuous  through  the  doublet  surface.  If  one  evaluates  the  normal 
derivative  at  the  centroid  of  the  doublet  surface  one  finds  that  the 
magnitude  decreases  to  zero  in  the  limit  as  this  finite  area  increases 
to  include  the  entire  plane. 
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2.T  Preasure  Coefficient  Evaluation 


Using  the  equations  Just  derived  for  the  influence  coefficiento 
it  is  poaaible  given  the  geometry  to  define  the  necessary  influence 
coefficients  needed  for  the  set  of  simultaneous  equations  which  con- 
sti^’ute  the  surface  boundary  conditions.  Once  the  solution  in  terms  of 
the  doublet  strength  la  known  the  local  velocities  may  be  determined 
essentially  by  back  substitution.  Given  the  velocities  we  may  then 
determine  the  local  pressures  using  the  equation  of  motion  given  by  equation 
(2.1.9),  that  is 


dt 


P 


(2.7.1) 


After  expanding  this  equation  using  vector  identities  and  interchanging 
vecxor  operations  we  find  that  the  equation  of  motion  may  be  written  as 


(2.7.2) 


Since  this  relation  must  hold  throughout  the  flow  field  it  is  necessary 
then  that 

y  ^  ^  ^  7'  •  /T'  _  (2.7.3) 


If  we  now  transform  this  equation  to  the  non-inertial  primed  reference 


system  we  find  that 
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j  -  /  -  / 


(2.T.>*) 


_  / 

Evaluation  of  the  constant  //  for  ^  allows  us  to  write 


■j 


-  /  -  / 


(2.7.5) 


where 


--  '  £^^  X  a'  ^  . 


(2.7.6) 


If  we  define  a  pressure  coefficient  such  that 


-  Z’ 

_i_ _ ijst _ 

J-  /<-'  ^ 


(2.7.7) 


we  find  then  using  equation  (2.7*5)  that 


-- 


-  / 

i/ 


(2.7.8) 


'hus  given  the  velocities,  the  pressure  coefficients  may  be  determined. 
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CHAPTER  3.  DISCUSSION  OF  THE  SOLUTION  SCHEME 

3.1  Method  of  Preaentation 

In  the  diacuaaion  to  follow  in  this  section  we  will  discuss  more 
specifically  the  actual  solution  scheme  employed  beginning  with  the 
description  of  the  rotor  blades  and  carrying  this  through  to  the 
determination  of  rotor  forces.  The  symbols  uced  will  in  most  cases 
correspond  to  those  employed  in  the  computer  program.  All  lengths  and 
velocities  are  taken  to  be  nondimensionalized  on  the  basis  of  the  rotor 
radius  (R)  and  tip  speed  (f*®)  unless  specified  otherwise. 

3.2  Rotor  Planform  Description 

We  have  chosen  to  model  the  rotor  geometry  in  terras  of  a  blade 
fixed  orthogonal  coordinate  system.  In  particular  we  have  assumed 
that  there  are  I  identictl  rotor  blades  whose  span  axes  are 

displaced  from  each  other  in  the  ^ y'J  plane  of  rotation  by  an 
angular  measure,  DT,  given  by 

''  .  (3.2.1) 

We  further  assume  that  the  blades  are  numbered  consecutively  in  the 
direction  of  rotation  such  that  blade  one  has  its  span  axis  projection 
coincident  with  the  y'  axis.  It  follows  from  this  that  only  one 
blade  need  be  defined  as  the  other  blades  may  be  defined  from  the 
geometry  of  blade  one  by  a  simple  axis  rotation.  Figure  2  visualizes 
the  coordinate  system  convention  indicating  also  the  relation  between  the 
blade  fixed  coordinate  systems  and  a  fixed  inertial  reference  system 


(' X,  ^  description  of  blade  one  we  have  allowed  for 

two  approaches.  The  first  approach  to  describing  the  rotor  blade 
geometry  is  a  completely  genereii  method.  In  essence  wc  assume  the 
rotor  blade  is  positioned  in  space  at  some  lifting  attitude.  We 
require  only  that  it  have  some  non  zero  hub  radius  and  that  its  span 
eucis  projects  along  the  y  '  reference  euls.  We  assume  thrn  that  the 
blade  can  be  described  by  a  set  of  grid  point  coordinates  in  the 

'  y  '  £  'J  system.  If  we  let  represent  the  number  of  span 

stations  and  and  ^  represent  the  niunber  of  upper 

sxirface  and  lower  surface  chord  stations  per  spfn  station,  then  the  blade 
one  description  involves  a  table  of  coordinates  such  that  we  specify 


for 

/V 

/V 

X 1  L 

--  I  , 

and 

/V 

on  the  upper  surface 

on  the  lower  surface 

where 

/'O 

-  1 

is  the  root  chord  section 

'  X7 

is  the  tip  chord  section 

// 

L 

is  the  leading  edge  span  section 

/V 

is  the  trailing  edge  span  section 

on  the  upper  or  lower  surface. 

We  will  use  XJ  d  etc.  as  distinct  from  etc.  to  clearly 


distinguish  the  blade  one  surface  grid  points  from  any  other  general 
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point.  Note  that  the  K  system  is  coincident  with  the 

primed  reference  />  '  y^'  /'J  system  as  shown  in  r’ijpire  3. 

A  second  method  for  describing  the  rotor  planfonu  which  we  have 
chosen  to  include  in  the  programmed  scheme  involves  the  description 
of  the  grid  points  for  a  linearly  twisted,  constant  chord,  constant 
airfoil  section  type  rotor  blade. 

In  this  method  we  assume  the  airfoil  section  is  a  NACA  OOXX 
series  airfoil  where  XX  represents  the  percent  thickness.  In  des¬ 
cribing  the  blade  we  now  specify  our  grid  points  as  follows : 

/V 

/V  /f  ^ ^  ^ 

for  /^  -  Ij  and  /V or  ^  where  SG 

represents  a  span  station  expressed  in  terms  of  percent  true  span  and 
CGU  or  CGL  represents  a  chord  station  expressed  in  terms  of  percent 
chord  on  the  upper  or  lower  surface  respectively.  The  remaining  para¬ 
meters  necessary  to  define  the  rotor  blade  surface  include  the  hub 
radius  (Bl),  the  chord  length  (C)  and  the  airfoil  percent  thickness 
ratio  (TC).  In  applying  the  linoar  twist  to  the  rotor  blade  the  geometric 
pitch  (theta)  at  any  span  location  is  given  by 

--  yy-  (3.2.2) 

where  TO  -  geometric  pitch  at  the  root  chord 

TT  -  amount  of  twist 


SG 


true  percent  span  station. 
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This  twist  is  applied  about  the  span  axis  which  we  assume  to  bv* 
located  behind  the  leading  edge  intersecting  the  chord  line  at  a  chord 
station  given  by  B3.  The  true  span  {B2)  is  taken  to  be  the  actual 
blade  span  length  and  is  given  by  the  difference  between  the  rotor 
radius  (R)  and  the  hub  radius.  Figure  3  attempts  to  visualize  these 
geometry  parameters. 

Given  this  information  we  can  now  calculate  the  grid  point  location 
in  terras  of  our  /V/,  YJ,  reference  system  and  express  the  results 
as 


yi 

/V 

^  1 

/y 

yi 

//  1^,-1 

which  results  in  a  surface  grid  point  description  compatible  to  the 
general  description  method  described  earlier.  We  hp./e  incorporated  in 
this  scheme  as  an  option,  a  method  to  fair  the  tips  by  simply  forcing 


the  surface  coordinate  at  the  tip  to  lie  along  the  chord  line. 
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3.3  Elemental  Planar  Surface  Description 

In  the  description  of  the  rotor  blade  surfaces  and  rotor  blade 
wake  surfaces  we  are  given  surface  grid  point  locations.  The  surface 
defined  within  four  control  points  is  a  warped  surface  ideally  des¬ 
cribing  the  blade  or  wake  surface.  Because  of  the  numerical  approxima¬ 
tions  applied  to  the  integral  boundary  condition  it  is  necessary  to 
define  a  planar  surface  as  our  control  surface  element.  In  a  manner 
simileur  to  that  found  in  Reference  lU  we  can  define  a  pleuiar  element  given 
the  four  grid  point  locations. 

We  assume  that  the  four  grid  points  whose  given  coordinates  are 
y'i),  are  numbered  clockwise  when  looking  down  on  the 
control  surface  from  the  external  flow  field.  We  can  now  construct 
two  tangent  plane  vectors  and  7^  such  that  their  cross 

product  will  define  a  unit  elemental  surface  normal  direction.  That 
is  let 


where 


/ 

-y- 

_  / 

a 

II 

■<J  ~ 

✓ 

■ 

- 

y; 

- 

4'- 

■'jy  -- 

- 

-  4 

(3.3.1) 
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then  /V  "  -5  ^  whose  components  sure 


A'y  -■ 

-  5'^. 

-- 

Ll  X 

is  Just 

// 

(3.3.3) 


>t  - 


//// 


(3.3.1*) 

The  elemental  surface  plane  can  be  specified  using  the  known  unit 


normal  and  a  calculated  average  point  ^ ^ to  be  in  the 


plane  auid  calculated  as  the  simple  average  of  the  corner  point 
coordinates.  In  order  now  to  project  the  given  grid  points  into  this 
elemental  plane  along  the  unit  normal  we  first  find  the  distsuice  J 


;th 


that  the  i  grid  point  is  removed  from  the  average  point.  This  is 


given  by 


(3.3.5) 


Then  the  projection  of  the  grid  point  into  the  plane  along  the  unit 


normal  results  in  a  set  of  grid  projected  coordinates  y' /  , 


given  by 


X/’,  - 

X/  / 

yy' 

y;, 

/y 

(3.3.6) 
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We  can  now  construct  an  elemental  fixed  axis  system  //  y,  T J  such 

that  the  X  axis  is  directed  along  the  unit  normal  /aV  .  the 
/  SLxis  is  directed  along  the  line  Joining  projected  corner  point 
1  with  projected  comer  point  U  and  the  ajcis  direction  is  given 

AS 

by  X  ^  . 

Let  vector  be  given  by 


y 


A  / 


/ 


(3.3.T) 


where 


(3.3.8) 

< 

-  Zl'  ^ 

Then  the  unit  vector  is  given  by 


/^  / 


--  y  ^  f 


f 


/ 


(3.3.9) 


Let  be  a  unit  vector  along  the  axis ,  then 


where 


(3.3.10) 
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ll 
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(3.3.11) 
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A 

If  we  now  take  the  origin  of  the  elemental  rj  coordinate 


ayatem  to  be  at  projected  corner  point  (l)  then  the  tremsformation 
equation  relating  the  body  fixed  coordinate  system  to  the  elemental 
coordinate  system  is  given  by 


JJ  ' 

X  * 

-  X/ 


(3.3.12) 

During  the  influence  calculations  which  will  be  described  later 
it  will  be  necessary  to  designate  one  point  on  this  elemental  control 
surface  at  which  the  influence  is  to  be  ceLLculated.  We  shall  take  as 
this  control  point  the  centroid  of  the  elemental  surface.  li  can  be 
shown  that  the  centroid  location  given  by 


where 


^4 

i 

--  -  i/D 


(3.3.13) 


(3.3.li4) 
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3.4  Blade  Influence  Calculations 

We  have  to  this  point  in  the  discussion  essentially  described 
the  geometry  of  blade  one  as  a  set  of  upper  and  lower  surface  grid 
points.  Furthermore,  we  have  presented  a  method  to  specify  an 
elemental  planar  control  surface  in  terms  of  its  location  in  apace 
with  respect  to  a  blade  fixed  reference  system  and  with  respect  to 
its  own  elemental  coordinate  system.  Therefore,  we  have  the  informa¬ 
tion  necessary  to  compute  the  influence  coefficient  of  one  blade  element 
on  smother  blade  element. 

In  order  to  make  the  problem  more  easily  tractable  we  shall  first 
define  a  new  indexing  system  for  identifying  the  particular  elemental 
surface.  Recall  that  in  the  description  of  blade  one  we  have  specified 
spam  stations,  and  ^  upper  and  lower 

chord  stations  where  grid  point  M^l,  N^l  corresponds  to  the  root  chord 
at  the  leading  edge  of  the  upper  or  lower  surface ,  and  grid  point 
yV  /y  ,  yV-'  yVL,  ^  corresponds  to  the  tip  chord  at  the 

trailing  edge  on  the  upper  surface.  If  we  index  the  elemeital  control 
surface  with  the  indices  taken  to  be  the  index  of  the  most  immediate 
inboard  span  station  and  forward  chord  station  then  essentially  we  can 
describe  the  particular  elemental  control  surface  as  (element If 
we  define 


-  Z 

(3. 1^.1) 

-  z 

(3.1<.2) 

-  / 

(3.»<.3) 
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then  we  ceui  state  that  there  are  elemental  control 

surfaces  on  the  upper  blade  one  surface  and  likewise,  there  are 

elemental  control  surfaces  on  the  lower  blade  one  surface. 
Because  we  will  be  specifying  the  influence  of  one  element  on  another 
element  it  behooves  us  to  index  the  elements  with  a  single  index.  We 
will  Identify  the  elements  as  where  the  elements  are 

numbered  consecutively  beginning  on  the  upper  surface  at  the  root 
leading  edge  element  and  proceeding  first  chordwlae  and  then  spanwise. 
If  J ■'  is  the  last  element  on  the  upper  surface  then 

j .  T  ^  /  is  taken  to  be  the  root  leading  edge  lower  surface 
element  auid  the  indexing  proceeds  to  in  a  similar  manner  to 

the  upper  surface.  It  follows  then  that  for  a  particular  control 
element  specified  as  //%/^J  then  in  terms  of  index  I  we  have, 


on  the  upper  surface  (3.1*.^4) 


^ J /  /V ^  on  the  lower  surface  (3.1*.5) 


where 


(3.14.6) 


The  computation  of  the  blade  influence  coefficients  now  proceeds 
in  a  straight  forward  fashion.  That  is  suppose  we  first  define  once 
and  for  all  the  elemental  control  surfaces  of  blade  one. 

Beginning  on  the  upper  surface,  say  elemental  surface  (M.H),  we  can 
immediately  using  equation  (3.1*.1»)  index  this  element  as  elemental 


% 


surface  I.  Now  proceeding  in  a  clockwise  fashion  around  the  element 
we  define  the  comer  point  coordinates  as 


/ 

•</  - 

il  ' 

(3.1*. 7) 
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XI 
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- 
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-  a  AO  . 

Note  for  the  lower  surface  element  we  would  define  the  corner  point 
coordinates  in  a  clockwise  fashion  as 


xy  - 

XJ  ro 

//  '- 

y  1 1  A^^A^ 

-- 

?  j  ^A'  ^ 

-- 

X  i  ^  /O 

//  - 
/ 

^a^ai,  /y 

XJ  ^  aaaS,  /'O 

(3.1*. 12) 


146 
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^  = 

- 

V  - 

//  - 

yj  /'j 

- 

^  ^  j  . 

(3.14.13) 


(3.14.114) 


Nov  using  the  equations  developed  in  the  discussion  of  the  elemental 
planar  surface  (section  3.3)  we  can  define  all  the  geometric  relations 
needed  for  later  influence  calculations.  That  is,  ve  develop  and  save 
the  following  parameters  for  all  elements  I  of  blade  one: 

A)  Transformation  matrix  relating  the  primed  reference  system 
to  the  elemental  reference  system. 

B)  Origin  of  the  elemental  coordinate  system  with  respect  to  the 
primed  reference  system. 

C)  Elemental  corner  point  coordinates  with  respect  to  the 
elemental  reference  system. 

D)  Elemental  centroid  coordinates  with  respect  to  the  elemental 
reference  system  and  with  respect  to  the  primed  reference 
system. 

E)  Elemental  area. 

Given  now  thr/  element  I  of  blade  1  on  which  we  wish  to  find  the 
influence  of  element  J  of  blewie  1  we  would  proceed  as  follows: 

a)  Recall  the  location  of  the  centroid  of  element  I  with  respect 
to  the  y  I  reference  system,  ie  y;  ^ 


0 


hi 


B)  Using  the  transformation  matrix  relating  the  primed  reference 

system  to  the  Jth  elemental  Xj^  coordinate  system, 

transform  the  location  of  the  coordinates  /'y I  //  ^ 

coordinates  / X, 

C)  Given  nov  the  corner  point  coordinateiof  element  J  with  respect 

to  its  coordinate  system,  ie  X/y/)  *  calculate 

using  equation  (2.6 .18) the  velocity  Influence  coefficients  in 
the  /y  ^  and  X  direction  at  point 

that  la,  define  ^  ^  and  _Xy,  ^ 

D)  Since  these  are  the  influence  coefficients  expressed  with 
respect  to  the  Jth  elemental  coordinate  system,  transform  these 
coefficients  to  the  Xx'y'^  7'J  system  using  the  inverse 

of  the  transformation  matrix  of  B  above  such  that  the  final 
result  is  in  the  form  of  velocity  influence  coefficients 
expressed  with  respect  to  the  primed  reference  system,  that 


is.  ^ 


Since  the  boundary  condition  requires  us  to  define  the  influence  co¬ 
efficients  on  all  X blade  one  elements  of  all  elements  of  all 
blades  ^ '  ^y  would  simply  repeat  the  procedure  outlined 

above  for  all  J  elements  of  blade  one,  and  for  exactness,  ve  should 
replace  the  Jth  subscript  above  with  J^  to  indicate  that  this  represents 
the  influence  of  the  Jth  element  of  blade  one  on  the  Ith  element  of 
blade  one. 

In  order  now  to  compute  the  elemental  influence  coefficients  of 
blades  /■  X ^  elements  I  of  blade  one  we  have  to  alter  the 


,  and 


procedure  outlined  earlier.  Note  first  that  we  have  asaxuned  all  rotor 
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blades  are  identically  described.  Thus  the  J  elements 1  surface  of 
blade  //I  similar  to  the  J  elemental  surfaces  of  blade  2  , 

Thus,  since  ve  are  given  the  geometric  parameters  defining  blade  one 
elements ,  all  that  ve  need  to  do  In  order  to  define  the  parameters  of 
the  J  elements  of  blade  2  /  J  la  to  apply  a  simple  transformation 
Involving  a  rotation  of  axis  about  the  2  or  21  axis.  The 
rotation  angle  la  of  course  the  angle  between  the  blade  span 

axis  and  the  blade  Z 2 2  span  eucis.  Once  this  transformation  Is 

accomplished  the  outlined  procedxrre  above  may  be  followed.  The  final 
res'olt  Is  then  a  set  of  Influence  coefficients  relating  the  Influence 
of  blade  L,  elemental  surface  J  on  the  bleide  one  Ith  elemental  surface 
centroid.  The  total  Influence  then  of  all  blade  J  elements  on  all  I 
elements  can  be  written  with  respect  to  the  primed  reference  system  as 


2-2  ^2,2^ 

(3.1*. 15) 

(3.1*.i6) 

(3.1*. 17) 


Thus  the  blade  Influence  coefficients  are  determined. 
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3.5  Wake  Influence  Coefficients 

The  problem  concerned  with  calculating  the  wake  influence  co¬ 
efficients  once  the  wake  geometry  is  established  is  essentially  the 
same  problem  encountered  with  calculating  the  blade  influence  co¬ 
efficients.  For  the  wake  influence,  we  are  concerned  with  calculating 
the  influence  on  blade  one  element  I,  of  a  wake  element  as  it  leaves 
the  trailing  edge  and  spirals  down  in  a  streanwise  direction.  Since 
we  have  shown  that  the  doublet  strength  remains  constant  in  the  wake 
along  a  streamline,  then  ve  may  numerically  integrate  the  influence  of 
streamline  segments  and  sum  the  results  to  present  the  net  effect  of  a 
particular  wake  streamline  on  an  element  I. 

In  this  problem,  because  of  the  numerical  approach,  we  are  concerned 
not  with  a  single  weike  streamline  but  rathei  k  '  a  series  of  adjacent 
streeiffilines  which  leave  the  rotor  blade  in  finite  segments  on  each  of 
which  the  doublet  strength  is  a  constemt  given  by  the  difference  in 
upper  euid  lower  blade  surface  trailing  edge  element  doublet  strengths. 

As  such  we  wish  now  to  calculate  the  influence  of  spiraling  finite  width 
wake  segments  on  blade  elements  1=1,  Imax* 

We  will  choose  to  represent  the  wfdie  as  I4M  spiral  segments ,  each 
of  which  trails  from  one  of  the  J  blade  trailing  edge  surface  elements 
previously  discussed.  In  order  to  demonstrate  the  feasibility  of  a 
lifting  rotor  surface  analysis  and  to  begin  the  free  wake  analysis  we 
will  choose  to  model  the  wake  in  terms  of  a  prescribed  classic  wake 
model.  The  description  of  this  wake  is  relatively  simple  and  is  used 
here  because  of  that  fact.  In  no  way  is  this  theory  limited  to  such  a 


wake. 
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Frwn  the  rotor  blade  one  descr  iption  we  know  the  tVax  trailing 
edge  span  station  coordinates.  These  may  be  expressed  in  terms  of  a 
primed  cylindrical  coordinate  system  such  that  the 

dimensional  coordinates  are  given  by 


j  /y  --I,  , 

(3.5.1/ 

It  follows  then  that  the  streamline  leaving  the  trailing  edge  at  euiy 
of  these  span  stations  will  follow  a  helical  path,  the  equation  of  which 

j  / 

is  expressed  in  terms  of  a  parameter  as  follows , 


-  constant 


where  is  a  constant  rotational,  rate  taken  as  positive  in  the  usual 

sense  about  the  p'  axis,  and  is  a  constant  velocity  taken  as 

positive  along  the  £  '  axis.  If  the  distances  are  nondlmenslonalized 
in  the  usual  sense  using  the  blade  radius  (P),  and  the  velocities  are  non- 
dlmensionalized  using  the  rotor  ratlonail  tip  speed  ^ ,  then 
the  above  set  of  equations  becomes 


where 


-  constant 

- 

(3.5.3) 

(3.5.**) 
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In  this  ancLlysls  we  will  take  the  velocity  ^2^' J  to  be  the  sum  of 
the  rotor  climbing  speed  and  a  specified  constant  increment  downwaish 
such  that  the  sum  represents  an  average  constant  momentum  downwash 
expressed  with  respect  to  the  blade  fixed  reference  system. 

Having  established  these  equations  it  is  now  a  simple  matter  to 

describe  the  wake  in  terms  of  elemental  surfaces  Just  as  was  done  for 

the  blade  surface.  We  will  define  the  wake  as  consisting  of  MM  spiral 

segments  such  that  the  Mth  spiral  wake  segment  has  as  its  origin  the 

trailing  edge  segment  lying  between  trailing  edge  grid  points  M  and 

M4’l.  Each  spiral  wake  segment's  influence  then  on  blade  one  element  I 

will  be  calculated  as  the  sum  of  a  finite  number  of  trailing  segments 

derived  by  incrementing  the  wedce's  displacement  from  the  blade  in  terms 

of  the  parameter  The  same  procedure  as  outlined  in  the  bletde 

influence  coefficients  section  3.^  is  followed  with  the  result  being 

expressed  as  a  set  of  Influence  coefficients  expressed  in  terms  of  the 

primed  ^ ^  U  reference  system  as  j  j  ^ 

'y  y  ' ■' 

where  7]/'  la  bhe  influence  coefficient  in  the  direction, 

derived  from  the  influence  on  blade  one  element  1  of  the  wake  M  segments 
trailing  the  ^ number  of  blades  from  the  Mth  trailing  edge 
elemental  surface. 

In  the  actual  computation  scheme  the  wake  influence  calculations 
are  terminated  depending  on  one  of  the  following  constraints: 

A)  After  a  finite  number  of  wake  turns  are  completed. 

B)  After  the  influence  coefficient  on  a  specific  blade  element 
I  has  reached  a  value  less  than  a  specified  percenteige 


amount  of  a  close  blade  wake  element. 
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C)  After  a  specified  absolute  number  of  vake  Influence 
*  calculations  are  performed. 

The  size  of  the  wake  elements  are  determined  by  the  length  of  the 
trailing  edge  span  stat:*^n  grid  points  and  the  angular  measure  , 

An  option  is  Included  in  the  program  to  allow  '  to  be  specified 

as  a  certain  value  until  a  given  number  of  wake  turns  are  completed  at 
which  time  will  take  on  a  second  specified  value.  This  allows 

for  finer  wake  elemented,  control  surfaces  to  be  specified  in  the  wake 
region  near  blade  one  and  coarser  wake  elemental  control  surfaces  to  be 
specified  in  the  far  wake  region. 
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3.6  Generation  of  the  Set  of  Surface  Boundary  Conditions  and  the 

Numerical  Solution 

To  this  point  in  the  discussion  the  determination  of  the  geometry 
of  the  blade  and  wake  system  and  the  analysis  of  the  entire  set  of 
influence  coefficients  have  been  presented.  We  are  concerned  now  with 
the  actual  formulation  of  the  set  of  surface  boundary  conditions.  Using 
the  results  of  the  previous  two  sections  we  may  write  the  Induced  normad 
velocity  on  the  centroid  of  a  blade  element  I  se 

-  /  _  jUj.  ^ 

y^j  ■  ^  4.JI  (3.6.1) 

where,  except  for  the  doublet  strength,  everything  else  is  analytically 
determined.  The  first  summation  term  represents  the  blade  Induced 
velocities  and  the  second  summation  term  represents  the  wake  Induced 
velocities.  The  element  I  unit  normal  expressed  with 

respect  to  the  primed  coordinate  system  has  been  previously  determined. 

flf  represents  the  normal  influence  coefficient  obtained  by 
vectorlally  dotting  the  respective  vector  primed  influence  coefficient 
with  the  unit  normad  Hj  . 

Once  the  doublet  strengths  have  been  cadculated,  the  surfau:e 
velocity  i>otentlad  is  computed  using  equation  {2.U.6).  The  induced 
surface  velocity  is  then  computed  by 

(3.6.2) 

as  / 

where  represents  the  tangential  surface  derivative. 


It  vu  shown  earlier  that  the  primed  free  stream  velocity  on  the 


centroid  of  an  elestent  I  can  be  vrittw  as 


(3.6.3) 


where  ^  is  the  rotor  rotational  rate  i  ')  ^  Vy  is  the  rotor 

climb  velocity  y' ^  is  the  primed  vector  to  the 

centroid  of  element  I.  Nov  ^  and  are  known, 

—  / 

thus  is  determined  and  it  follows  then  that  the  primed  free 

stream  velocity  normal  to  the  surface  element  I  given  by 


(3.6.U) 


is  determined. 

Substituting  these  results  into  equation  (3.6.2)  it  follows  then 
that  the  set  of  surface  boundary  conditions  to  be  satisfied 
is  given  by 


Since  was  shown  in  section  3.5  to  be  related  to  specific 

values  the  set  of  equations  above  are  in  fact  in  the  form  generally  given 

by 


l^J  Jx/  =  /£/ 


(3.6.6) 


55 


wh«re  the  matrix  A  la  an  square  matrix,  X  is  a 

column  matrix  and  B  la  an  column  matrix.  We 

have  identified 

✓  -  - 

^7  jfTJ  (3.6.7) 

✓ 

^  (3.6.8) 

The  solution  to  the  simultaneous  set  of  linear  algebraic 

equations  la  carried  out  numerically  by  either  of  two  methods  depending 
on  the  program  option  selected.  One  method  employs  the  Gaussieui  eli¬ 
mination  scheme  with  pivoting.  This  method,  which  la  discussed 
completely  in  Reference  17,  is  an  exact  solution  method  whose  scheme  of 
operation  is  to  successively  convert  the  defined  matrix  into  a  matrix  of 
one  less  row  euid  column  until  in  effect  an  upper  trlanguleur  matrix 
results.  The  solution  is  obtained  then  by  a  back  substitution  scheme. 
Interchanging  of  rows  is  done  in  order  to  make  the  pivot  diagone^  terms 
dominate  idilch  is  a  condition  necessary  to  ensure  accuracy. 

A  second  method  chosei  to  obtain  a  solution  to  the  set  of  equations 
is  aui  iterative  approximate  method  known  eis  "The  Method  of  Simultaneous 
Displacement  for  Linear  Systems".  A  complete  discussion  of  this  method 
is  presented  in  Reference  18.  In  essence,  however,  this  scheme  involves 
expressing  the  linear  equation  in  the  form 


/X/  ^  /^7/x/  /  /X/ 


(3.6.9) 


>0 


vhtre  B  iB  an  aquare  matrix  and  C  Is  the  given  column 

matrix.  An  initial  approximation  vector  is  taken  aa 


Xr  - 


(3.6.10) 


Successive  approximations  Xj  are  generated  by  the  Iteration 

J 


(3.6.11) 


This  continues  until  the  criterion 


1  ^  / 
/  Xj  -  Xj  / 

/  x/^^  / 


X  X  X^A  X~ 

(3.6.12) 


is  satisfied. 

In  actual  practice  the  iterative  scheme  is  the  solution  of  choice 
for  large  sets  of  simultaneous  equations.  This  point  will  be  discussed 
In  a  later  discussion  section. 
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3.7  Detenalnation  of  the  Velocities,  Pressures,  and  Forces 

Once  the  solution  in  terms  of  the  doublet  strengths  are  known 
it  is  a  simple  matter  to  determine  the  Induced  primed  velocity  on 
all  blade  one  elements  using  equation  (3.6.2).  The  pressure  coefficient 
at  the  centroid  of  element  I  is  then  given  by 


X 


(3.7.1) 


where  ^  is  taken  here  to  be  the  ratio  of  the  rotor  climb  speed  to 
the  tip  speed  r)  . 

The  non-dimensioneLl  force  acting  on  any  element  I  is  given  by 


S' 


(3.7.2) 


v^/  ^  ^  f 

Where  S  la  any  y  or  ^  direction. 

The  non-dimensional  plane  of  rotation  torque  acting  on  any  element 
I  is  given  by 


(3.7.3) 


The  total  rotor  axial  thrust  coefficient  is  given  by 


(3.7.14) 


M 
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Similarly  the  total  rotor  torque  coefficient  is  given 


2  77 


(3.7. 


5) 


This  concludes  the  step  by  step  discussion  of  the  solution  scheme. 
As  a  way  of  a  short  summary  at  this  point  we  have  shown  how  if  we  are 
given 

A)  The  rotor  geometry  including  specific  surface  grid  points. 

B)  The  rotor  non-dimensional  climb  ratio  ^  A  “  ^ 

C)  Some  prescribed  wake  geometry  specifically  taken  in  the 
previous  discussion  to  be  a  simple  classical  helix  shape 

then  we  can  calculate  the  pressure,  velocities,  and  forces  acting  on 
the  blade. 


A 
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CHAPTER  1*.  DISCUSSION  OF  THE  COMPUTER  PROGRAM 
U.l  General  Description 

In  order  to  demonstrate  the  feasibility  of  this  lifting  surface 
theory  and  its  application  to  the  prediction  of  rotor  blade  flow  field 
analysis  a  computer  program  was  written.  The  theory  as  described  in 
the  previous  sections  is  a  relatively  straightforward  theory,  however, 
its  application  becomes  rapidly  complicated  due  to  the  geometry  and 
large  system  of  elemental  control  surfaces  necessary  for  accurate  results. 

As  a  means  of  perceiving  the  extent  of  the  problem  we  may  consider 
at  this  point  a  reasonable  end  problem  that  we  might  wish  the  computer 
program  and  computer  to  handle. 

Suppose  we  are  given  a  two  bladed  rotor  system  whose  geometry  is 
specified  along  l6  spajti  stations,  21  upper  surface  chord  stations  and 
21  lower  surface  chord  stations.  This  implies  that  we  will  eventually 
describe  the  three  component  vector  velocities  on  600  basic  elemental 
surfaces.  We  have  seen  that  the  description  of  esuih  elemental  surface 
alone  requires  at  least  25  parameters  which  for  this  posed  problem  means 
a  total  ofl5|000  such  parameters.  The  blade  influence  calculations  alone 
total  720,000.  If  we  describe  the  wake  with  a  total  of  1155  basic  wake 
elements  we  require  a  total  of  1,386,000  wake  influence  calculations  to 
be  performed.  The  simultaneous  set  of  equations  to  be  solved  is  a  600  x  600 
system. 

The  computer  progreun  was  run  for  this  problem  and  will  accommodate 
up  to  2000  basic  bleuie  elements  which  is  the  critical  size  determining 
factor.  It  should  be  obvious  that  the  constraints  to  this  type  of 
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computer  solution  vlll  be  the  available  size  of  the  main  core  computer 
storage  and  the  computatlonaj.  speed  of  the  computer  itself. 

The  program  was  designed  to  be  run  on  West  Virginia  University's 
IBM  360  Model  75  computer.  The  greatest  amount  of  computer  core 
stbrage  required  at  any  one  time  is  21*0,000  bytes  with  U  lytes  being 
required  for  a  word  length.  The  program  requires  a  minimum  of  three 
accessible  external  storage  devices  with  at  least  two  of  these  devices 
being  a  high  speed  storage  device.  At  WVU  this  requirement  was  met  by 
using  the  system  disc  files  with  two  standeurd  length  nine  track  tapes 
being  used  for  permauient  storage .  Seven  track  tape  devices  sure  in¬ 
compatible  with  the  progrsun  as  written.  The  progreun  Itself  is  written 
in  FORTBAN  IV  langusige  and  was  compiled  on  the  IBM  FORTRAN  G  version 
compiler.  The  program  is  a  research  orientated  program  and  does  not 
necessaurlly  reflect  the  most  optimum  design  orientated  program,  that  is, 
it  contains  various  checks  emd  options  which  sure  not  absolutely  necesssury 
to  the  progrsun  solution.  We  recognize  also  that  there  exist  sureas  in 
the  program  in  which  the  exswrt  analytic  expressions  may  be  approximated 
and  advantage  of  these  time  saving  methods  have  been  used  in  this 
program.  However,  it  is  the  intent  of  this  program  not  to  Incorporate 
all  these  features  at  this  time.  The  philosophy  under  which  we  have 
written  this  program  is  to  msike  it  as  exact  a  numerical  scheme  as 
possible  so  that  if  an  infinite  number  of  control  elements  were  tsdcen 
the  solution  would  be  as  exaurt  as  possible  within  the  confines  of  any 
inconq)ressible  potential  flow  analysis. 

The  program  is  designed  on  a  modular  basis  with  each  module  being 
designed  on  the  premise  that  given  certain  input  informtlon  it  is 
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the  purpose  of  that  module  then  to  specifically  calculate  a  certedn 
aspect  of  the  overall  problem  and  to  supply  certain  output  information. 
The  advantage  here  is  that  it  allows  for  easy  internal  modification 
by  simply  replacing  this  module  with  another.  The  actual  sharing  of 
information  is  done  for  the  most  part  through  the  use  of  common  machine 
storage  centers.  In  this  way  the  amount  of  in  core  storage  necessary 
is  minimized  as  storage  locations  eu*e  used  at  different  times  for 
different  parameters. 

k.2  Specific  Program  Description 

The  overall  program  is  subdivided  into  seven  main  programs 
identified  as  follows: 

1.  tiAIN  -  INFLUENCE  COEFFICIENTS 

2.  MAIN  -  CREATE  SOLUTION  FILE 

3.  MAIN  -  SOLVE  BY  ELIMINATION 

U.  MAIN  -  SOLVE  BY  ITERATION 

5.  MAIN  -  TRANSFER  SOLUTION  FILE 

6.  MAIN  -  VELOCITY  CALCULATIONS 

7.  I-IAIN  -  PRINT 

The  purpose  of  each  program  is  generally  indicated  by  the  assigned  name. 
The  Print  program  is  actually  an  auxiliary  program  whose  sole  purpose 
is  to  retrieve  auiditional  non-essential  information  stored  on  the  two 
permanent  files.  It  was  used  during  the  program  check  procedures  to 
study  individual  elemental  influence  coefficients. 
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The  reason  a  series  of  main  programs  are  used  rather  than  one 
single  main  program  is  because  the  overall  program  logically  divides 
itself  into  these  areas  and  secondly  ve  can  minimize  the  amount  of 
external  storage  devices  and  internal  core  storage  required  dxiring 
any  one  series  of  computati«)nal  steps  thereby  allowing  for  the  computer 
to  be  run  in  a  time  sharing  mode  as  is  common  procedure  at  most  in- 
stedlatlons.  A  generad  computer  program  flow  diagram  is  given  in 
Figure  U  indicating  the  major  steps  in  the  program.  The  solid  lines 
indicate  continuation  to  the  next  step  vLtb  the  broken  lines  indicating 
continuation  into  another  main  program. 

In  Appendix  D  we  have  discussed  the  internal  structure  of  the  pro¬ 
gram  at  the  level  of  the  subroutine.  Basically  we  have  presented  each 
subroutine  as  an  entity  and  described  its  function  and  options.  Except 
for  the  MAIN-INFLUENCE  CALCULATIONS  program,  the  actual  program  operation 
is  straightforward.  We  have  presented  a  detailed  flow  diagram  for  the 
MAIN-INFLUENCE  CALCULATIONS  program  in  order  to  serve  as  an  aid  in 
discussing  this  program  operation  in  terms  of  its  various  options. 

We  have  also  included  in  Appendix  D  a  discussion  of  the  input  procedures 
and  output  information. 
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CHAPTER  5.  RESULTS  AND  CONCLUSIONS 

In  this  section  we  will  present  a  sujnmary  of  the  results  of  two 
major  case  studies.  These  cases  were  chosen  so  as  to  verify  the  mathe- 
matlced.  model  and  the  computer  program  for  the  case  of  a  prescribed 
wake  model.  Recall  that  the  present  work  is  not  directly  concerned 
with  the  force  free  wake  analysis  which,  in  fact,  would  require  subse¬ 
quent  wake  iterati  ’«  computer  rtins  to  be  made  for  each  of  the  cases 
presented.  The  results  presented  here  thus  represent  the  starting  or 
zeroth  order  solution  to  the  overall  rotor  analysis.  Subsequent 
solutions  will,  however,  make  use  of  the  same  program  with  the  prescribed 
wake  geometry  having  been  predicted  from  the  previous  program  solution. 

In  Table  1  we  have  presented  a  summary  of  the  rotor  description, 
elemental  blade  surface  description,  elemental  wake  surface  description 
and  free  stream  conditions  as  used  in  the  first  case  study.  This  single 
bladed  rotor  system  does  not  correspond  to  any  real  rotor  system  but 
was  chosen  because  its  analysis  requires  a  relatively  small  number  (l6o) 
of  basic  elemental  surfaces  for  each  blade  description.  Figure  5 
attempts  to  visualize  the  density  distribution  of  the  elemental  surfaces 
which  describe  the  rotor  blade  surface.  This  plot  is  a  planform  projec¬ 
tion  of  the  symmetrically  described  00  upper  and  00  lower  elemental 
surfaces.  We  have  spaced  the  elemental  surfaces  more  closely  in  those 
regions  where  we  expect  the  pressure  to  vary  most  rapidly.  The  rotor 
hub  radius  of  this  two  bladed  system  is  taken  to  be  90%  of  the  rotor 
radius,  thus,  the  rotor  blade  occupies  the  outboard  10^.  In  Figures 
6  through  9  we  have  presented  a  plot  of  the  pressure  coefficients  vs. 
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chord  station  at  each  of  the  tovor  span  segments  as  identified  on 
Figure  5.  The  pressure  coefficient  is  defined  on  the  basis  of  the 
rotor  rotational  tip  speed  squared  and  each  point  plotted  represents 
an  elemental  centroid  point  at  which  the  normal  surface  boundary  condi¬ 
tion  waia  satisfied.  These  plots  represent  a  realistic  chordwise  variation 
in  pressure  suid  show  the  three  dimensional  character  of  the  flow  near 
the  rotor  blcde  tips.  In  Table  2  we  have  summarized  the  computer  time 
required  for  this  case. 

In  order  to  more  clearly  show  the  three  dimensional  effects  we 
studied  the  same  two  bladed  rotor  system  eis  described  above  except 
that  the  geometric  pitch  angle  at  all  span  positions  was  now  taken  to 
be  zero  degrees.  All  other  conditions  remained  the  same  as  given  in 
Table  1.  This  einalyals  thus  corresponds  to  a  three  dimensional  non- 
lifting  potentled  flow  thickness  problem  and  as  such  no  wake  analysis 
was  performed.  Figure  10  is  a  plot  of  the  surface  velocity  vs.  chord 
station.  The  surface  velocity  here  plotted  is  the  velocity  at  each  of 
the  elemental  centroid  locations  of  the  tip  span  segment  (span  segment 
U  of  Figure  5)  divided  by  the  local  free  stream  velocity  given  by 
whore  Is  the  position  vector  to  the  elemental 

centroid  from  the  blade  center  of  rotation.  On  this  plot  is  also  shown 
two-dimensional,  zero  angle  of  attack,  NACA  0012  sectional,  potential 
flow  derived  surface  velocities  non-diraensionallzed  by  the  free  stream 
velocity  as  given  in  Reference  19*  Only  the  upper  surface  velocities 
are  plotted  since  the  velocities  are  the  same  on  the  lower  surface 
for  this  symmetric  airfoil.  Note  that  the  three  dimensional  blade 
velocities  are  less  in  magnitude  than  the  comparable  sectional  airfoil 
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vtlocitles  and  the  velocity  curve  defined  by  the  rotor  blade  points 
has  a  near  zero  slope  over  much  of  the  chord  distance.  These  differ¬ 
ences  represent  here  the  change  in  flov  character  between  a  three 
dimensioned  body  executing  a  constant  rotating  motion  compeu*ed  to  a 
two  dimensional  body  executing  a  constant  translatory  notion.  The  fact 
that  the  rotor  velocity  ratio  does  not  equal  the  two  dimensional  result 
is  due  to  the  three  dimensional  divergence  effect  on  the  flow  over  this 
small  blade  segment . 

The  second  major  case  study  is  an  analysis  of  an  actual  rotor  system 
for  which  experimental  surface  pressure  data  is  available  (reference  20). 
In  Table  3  we  have  presented  a  summary  of  the  rotor  description,  elemental 
blade  surface  description,  elemental  waJte  surface  description  and  free 
streaun  conditions  as  used  in  this  second  case  study.  The  density 
distribution  of  the  planfcrm  projected  elemental  surfaces  is  shown  in 
Figure  11.  Note  that  there  are  on  each  of  the  two  blades  600  basic 
blade  elemental  surfaces  symmetrically  distribv.ted  on  the  upper  and 
lo\rer  blade  surfaces  along  the  15  Identified  span  segments.  The  surface 
elements  are  more  concentrated  near  the  leewling  edge  where  we  expect 
the  surface  pressure  to  vary  most  rapidly.  In  Figures  12  through  22  we 
have  plotted  the  difference  in  the  upper  and  lower  surface  pressure 
coefficient  vs.  the  chord  station  at  11  selected  span  stations.  Kach 
point  plotted  represents  the  pressure  coefficient  difference  obtained 
at  two  respective  elemental  centroid  locutions  where  the  normal  surface 
boundary  condition  was  satisfied.  The  pressure  coefficient  is  defined 
here  on  the  basis  of  the  rotor  rotational  tip  speed  squared.  In  order 
to  better  appreciate  the  chordwise  and  spanwise  pressure  variation 
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we  have  presented  In  Figure  23  a  composite  curve  again  showing  the  pressure 
difference  vs.  chord  station  at  the  seune  selected  span  stations  In  dimen¬ 
sional  units.  Figure  2k  presents  an  experimentally  determined  plot  of 
the  absolute  upper  and  lower  surface  pressxire  difference  on  the  same 
rotor  as  a  function  of  chord  station  at  various  span  stations.  This 
figiire  Is  a  tracing  of  a  figure  presented  in  reference  20,  in  which  the 
pressure  difference  was  obtained  experimentally  at  6  chord  stations  at 
each  of  8  span  selected  stations.  Notice  that  the  same  scales  have  been 
used  for  both  curves  In  order  to  facilitate  comparison  between  the 
theoretical  and  experimental  results.  The  span  stations  do  not  correspond 
however,  so  care  should  be  taken  when  comparing  specific  values.  There 
Is  a  very  close  similarity  in  the  distribution  of  the  chordwlse  pressure 
between  the  experimentally  measured  and  computed  results.  Since  pressure 
taps  were  not  located  near  the  nose  of  the  experimental  rotor  the  authors 
of  reference  20  rounded  off  their  curves  rather  abruptly  In  these  regions 
which  probably  contributed  to  their  not  checking  as  closely  the  Integrated 
experimental  force  measurements  with  the  test  stand  value.  It  is 
Interesting  to  note  that  both  the  experimental  and  theoretical  results 
show  qualitatively  similar  losses  in  pressure  near  the  tip. 

The  Integrated  results  of  Figures  23  and  2k  are  shown  in  Figure  25. 
This  figure  shows  a  comparison  of  the  spanwise  loading  distribution  between 
the  experimental  results  and  the  theoretical  calculations.  The  agreement 
is  remarkably  close  in  view  of  the  fact  that  the  classical  prescribed 
wake  was  used  for  the  theoretical  lifting  surface  calculations.  As  shown 
in  reference  21,  the  prescribed  classic  wake  tends  to  predict  higher 
thrust  coefficients  than  the  experiment  shows.  This  reference  also 
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points  out  that  force  free  vokes  reduce  the  spamvise  loadln/^  over  the 
cluilc  voke  prediction.  The  present  results  are  therefore,  in  agreement 
and  it  is  anticipated  that  the  inclusion  of  a  force  free  vake  progreun 
will  bring  the  theory  into  better  agreement .  The  computed  thrust 
coefficient  for  this  case  was  Cr  *  O.OOUU  conipau*ed  to  the  experimental 
value  of  Cjm  0.0038.  Therefore  the  classic  rigid  wake  theoretical 
result  was  approximately  15t  high.  The  corresponding  theoretical  torque 
coefficient  was  computed  to  be  Ca»  0.00011  for  this  case. 

Figures  26  and  27  present  plots  of  the  predicted  relative  surface 
velocity  direction  on  the  upper  and  lower  rotor  surfaces  respectively 
for  case  2.  In  these  plots  we  have  forced  the  rotor  surface  to  conform 
to  the  plane  of  the  figvire.  The  direction  of  the  local  surface  velocities 
are  shown  drawn  os  a  constant  length  vector  from  a  series  of  surface 
control  points  which  approximate  selected  elemental  "it-.-itrold  locations. 

The  plots  ore  essentially  comparable  to  what  would  result  from  on 
experimental  flow  visualization  tuft  study.  The  figures  cl-srly  show 
the  three  dimensional  nature  of  the  flow  on  the  rotor  bl..  le  surface. 

Rote  the  circular  nature  of  the  flow  which  is  evident  at  all  span  stations 
but  which  la  most  prominent  at  Inboard  span  segments.  This  result  is, 
of  course,  expected  since  the  free  stream  component  is  truly  circular 

=  —  u)  R' )  .  However,  if  one  studies  these  results  more  clearly 
it  becomes  evident  that  the  Induced  velocity  component  does  alter  the 
flow  direction.  If  we  restrict  the  following  discussion  to  the  in¬ 
duced  spanwlse  surface  velocity  component  the  results  show  that  beyond 
approximately  the  50^  span  station  the  induced  velocity  is  directed 
Inboard.  Inboard  to  the  50/J  span  station  the  induced  velocity  la 
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directed  outboard  at  surface  points  forward  of  the  quarter  chord 
approximately  and  inboard  beyond  the  quarter  chord .  The  flow  deviation 
between  the  top  and  bottoa  surfaces  neeo*  the  tip  trailing  edge  is  of 
interest.  The  flow  leaves  the  top  surface  trailing  edge  at  approximately 
a  20  degree  difference  from  the  corresponding  flow  leaving  the  bottom 
surface.  The  directions  indicate  that  the  shed  vortex  sheet  is  already 
beginning  to  roll  up.  In  Table  U  we  have  summarized  the  canputer  time 
required  for  each  of  the  major  steps  in  the  execution  of  this  600  case 
study. 

Prom  the  summary  of  results  presented  we  conclude  that  the  lifting 
surface  theory  as  developed  is  applicable  to  rotor  system  analyses. 
Purthennore  the  application  of  this  theory  is  feasible.  The  results, 
however,  are  highly  dependent  on  the  prescribed  wake  model  used  in 
the  analysis  and  use  of  the  classical  wake  model  is  insufficient  to 
yield  accurate  results .  It  appears  that  the  prescribed  wake  geometry 
must  be  simlleu:  to  the  actual  free  wake  geometry  in  order  to  obtain 
accurate  loading  results. 

There  is  a  lack  of  experimental  data  or  theoretical  analyses  which 
would  help  us  confirm  the  theory  developed  here.  To  the  authors' 
knowledge  there  exist  very  little  available  literature  giving  rotor 
blade  surface  velocities  and  directions  and  loading  measurements  on  an 
''ting  rotor  system.  The  conclusions  we  have  drawn  above  have  as  a 
necessarily  assumed  ftie  present  theory  and  application  to  be 
The  general  potential  theory  certainly  has  been  proven  by  past 
application  on  planau*  body  motion  analyses  and  we  have  systematically 
checked  all  application  procedures  in  order  to  ensure  accuracy. 
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CHAPTER  6.  EXTENSIONS 

The  theory  and  formidable  program  presented  in  essv-nce  completes 
the  basic  approach  to  theoretically  analysing  a  lifting  rotor  using 
potential  lifting  surface  theory.  Much  work  needs  to  be  do.ie,  however, 
before  the  ajiplication  of  this  theory  can  be  used  routinely  for  rotor 
design  anedysis.  We  recognize  that  extensions  to  this  work  should 
proceed  along  two  major  routes.  One  major  route  of  study  should  be 
concerned  with  the  wake.  At  present  the  program  (not  the  theory) 
uses  a  classical  helix  wake  washed  avay  at  some  constant  velocity 
from  the  rotor  blades.  As  is  noted  by  various  authors  (explicitly 
cited  earlier  in  the  report)  the  wake  geometry  plays  a  major  role  in 
detenaining  blade  surface  velocities.  Obviously  a  more  accurate  pre** 
scribed  wake  representation  needs  to  be  incorporated.  At  the  present 
time  ve  are  extending  this  work  in  order  to  Incorporate  a  theoretical 
force  free  wake  emalysis.  This  extension  is  to  be  accomplished  by 
predicting  a  new  wake  geometry  configuration  based  upon  the  solution 
obtained  by  this  present  work.  In  this  manner  successive  solutions 
will  be  obtained  by  iteration  on  the  wake  geometry  such  that  the  force 
free  waUte  solution  will  result. 

A  second  major  area  where  further  extensions  to  this  work  is 
needed  is  in  the  area  of  reducing  the  required  computer  time.  This 
may  be  practically  accomplished  in  basicaLlly  two  ways,  that  is,  first 
reduce  the  total  number  of  influence  calculations  to  be  performed 
and,  second,  reduce  the  time  required  for  each  set  of  influence 
calculations.  The  former  may  best  be  accomplished  by  further  studying 


70 


the  Individual  wake  elemental  influence  and  optimizing,  based  on  this 
study,  the  actual  extent  and  elemental  description  of  the  wake  itself. 

The  latter  above  may  best  be  accomplished  by  using  approximate  relations 
to  describe  the  influence  coefficients  at  a  far  removed  surface  point 
instead  of  using  the  exact  relations  which  require  perhaps  twice  as  much 
con^juter  time.  The  numerical  results  shown  here  have  already  incorporated 
this  simplification  in  the  program. 

There  are  of  course  other  means  to  reduce  the  required  computer  time. 
Linearizing  the  surface  boundary  condition  would  reduce  the  time  required. 
Furthermore,  intorr>orating  an  experimentally  determined  wedce  geometry 
would  result  in  a  r-abstantlal  overall  computer  time  reduction.  This 
type  of  extension  would  make  the  present  work  suitable  for  studies  on 
existing  rotor  systems  but  not  necessarily  applicable  to  new  configuration 
rotor  systems. 

With  the  advances  b*lng  made  in  computer  technology,  especially  in 
the  areas  of  increasing  core  storage  coupled  with  the  extensions  above,  it 
does  appear  that  x..ic  rotor  lifting  surface  theory  will  allow  the  aerodynanicist 
to  perform  a  completely  theoretical  design  analysis  of  any  given  rotor  system. 

The  problem  which  remains  to  be  analysed  is  the  problem  posed  by 
a  helicopter  which  is  climbing  as  well  as  advancing  perpendicular  to 
its  axial  climb  direction.  The  analysis  of  this  rotor  system  requires 
an  unsteady  analysis  as  the  blade  loadings  now  become  functions  of  the 
blade  azimuth  position  thereby  increasing  the  number  of  unknowns  in  the 
problem.  The  theory  is  simply  an  extension  to  the  present  problem, 
however,  the  actual  solution  scheme  becomes  rather  complex. 
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AXIS  CONVENTION  RELATING  THE  POTENTIAL  AT 
A  FIELD  POINT  INDUCED  BY  A 
DOUBLET  SURFACE  DISTRIBUTION 
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ROTOR  SYSTEI-i  AXIS  CONVENTION 
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PLANE  OF  ROTATION 
FREE  STREAM  VELOCITY  PARALLEL  TO 
PLANE  OF  ROTATION 


FIGURE  2 


Illlllllll 


COMPUTER  PROGRAM  FLOW  DIAGRAM 


(Continued  from  previous  page) 
I 


I 


I 

I 


FIGURE  U  (b) 


ITHTICSN 


CHDRO 


FIGURE  12 


^HnRl>  STRTION 


91 


FIGURE  16 


CHaRI>  STRTIDN 


. 

J.  n 


ul 


I  I  I 

1 1  I ' 


I  lil 


n  j' 


u  ^ 


i  1  I  i  .  ;  i  i  i  1  •  i  i  i. _ i _ i-.-i 

a*/--  ffi-  9‘i-  i.’i-  z’l-  iri-  ff-  9- 


4  I  4  i 


cy3Mu“ij  d:^  -  ca3ddrd  d:; 


FIGURE  17 


FIGURE  l8 


FIGURE  20 


CMDRO  ^TRTI 


FIGimE  21 


97 


FIGURE  22 


RESULT  or 

VS.  CHORD  STATION 


STATION 


EXPERIMEinALLY  DETERMIMED 
PRESSURE  DIFFERENCE  VS  CHORD  STATION 
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UPF’EIR  RnTDR  5URF“RCE:  VE:i_nCITY  £>IRE:CT!aN  Rl_DT 


L-QWEIR  ROTOR  SURF'RCO  VEZLOCITY  OfRETCTIDN  Rl_OT 


l^OTOH  INPUT  ;:UMr^RY 
FOR 

CAUi-:  OTIDY  1 

HO'L'OR  DKCCHIPl'IOn 

NACA  OOi:'  AIRFOIL  SECTlOil 
HOTOK  RALiu:; 

ROTOR  IlUii  JIADIUO 
'IRUE  LPAW  LiENOTH 
COilUTAIlT  CHORD  LEiJGTll 
CPAiJ  AXIC  CJIORD  LTATIOH 
ROOT  CHORD  GEOilETRIC  PITCH 
LINEAR  TWICT  ABOUT  GI'AN  AXIL 
.iUi'IDER  OI’  liOTO.i  DIADEL 
.li'ITlIEi^  TIl'  [L.  FAIRED 

FREE  LTREAM  CONDITIONS 

HOVER  MODE 

ALLUMl^D  CONF.TAilT  V/AlJi  DOVAJWALH 

ELEMEiJTAL  BLADi;  fiURFACE  DELCRIPI’ION 

HUI4BF:ii  01'’  L.PAU  LIlGMiRTTL 
HUMHEH  Oh'  CHORD  LEGMENTL 

a)  UPPER  LIJRFACE  ONLY 

b)  LOWER  LURFACE  ONLY 
KDMENTAl,  AREAS  AJ(E  .'’.YiLMETRICAl.I.Y 

DELCHII3ED  WTTil  RE'U’iX’T  TO  CHORD  PLANE 

TOTAL  HUMBER  OF  liALIC  BLADE  ELEMENTS 
ELEMENTAL  WAKE  SURFACE  DESCRIPTION 


i.n 
0.9 
0.  1 
0 .  1 

10.0^ 

o.o'’ 
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V,.,^R  =  0 
wAa^R  =  0.05 
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)|0 
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PRESCRIBED  CLASSIC  WAKE  MODEL 

NUMBER  OF  SPAN  SEGMENTS  ^4 

NUMBER  OF  GTREAMWISE  SEGMENTS  85 

WAKE  ELElvlENT  ANGULAR  INCREI'lENT 

a)  NEAR  BLADE  REGION  20° 

b)  FOR  BLADE  REGION  (>80°  FROM  BIWYDE)  30° 

NUMBER  OF  WAKE  TURNS  COm’LETED  8.i 

TOT’AL  NUMBER  OF  BASIC  WA;CE  ELEI4ENTS  3*^0 


TABLE  1 
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COMPUTEH  KXECUTION  TIME  SUMMARY 


FOR 

CASi;  STUDY  1 

TIME 

EXECUTABU'J  SECTION  (SECONDS) 

I.  INFLUENCE  CALCULATIONS  I’HOdRAM 

A.  SLADE  INFLUENCE  CALCULATIONS  07 

(TOTAL  NUMBER  =  25,600) 

B.  WAKE  INFLUENCE  CALCULATIONS  228 

(TOTAL  NUMBER  =  ‘i*i,UOO) 

C.  O'riQ'JR  ROUTINES  6 

II.  CREATE  SOLUTION  FILE  PROGRAM  11 

III.  SOLVE  BY  ELIMINATION  PROGRAM  10** 

(EXACT  SOLUTION  METHOD) 

IV.  TRANSFER  SOLUTION  FILE  PROGRAM  1 

V.  VELOCITY  CAliCULATION  PROGRAM  1^ 

TO'rAL  ABOVE  '»5T 


TABLE  2 
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ROTOH  INPU'x'  SUMMARY 
FOR 

CAGE  STUDY  2 


ROTOR  DEGCRIPriQH 


IJACA  0015  AIRFOIL  GECTIOU 

ROTOR  RADIUS  1.0 

ROTOR  HUB  fiADIUS  0.1625 

TRUE  SPAR  LENGTH  0.ft3T5 

COii START  CHORD  LENGTH  0.1 

SPAI^  AXIS  CHORD  STATION  0.270B 

ROOT  CHORD  GEOMETRIC  PITCH  6.0° 

LINEAR  TWIST  ABOUT  SPAN  AXIS  0.0 

NUMBER  OF  ROTOR  BLADES  2 

NEITHER  TIP  IS  FAIRED 


FREE  STREAM  CONDITIONS 


liuVER  MODE  V.^ViLiR  =  0 

RE2.UMED  CONSTANT  WAKE  DOWNWAGH  w/<DH  =  Q.OUhU 

ELEMENTAL  SLADE  SURFACE  DESCRIPTION 

NUMBER  OF  SPAN  SIT/MENTS  15 

HUMBER  OF  CHORD  SEGMENTS 

a)  UPPER  SURFACE  ONLY  20 

b)  LOWIOR  SURFACE  ONLY  20 


ELl'MENTAL  AREAS  ARE  SYMMETRY TL! '' 

DESCRIBED  WITH  RESPECT  TO  CHORD  PLANE 

TOTAL  NUHIU'R  OF  BASIC  BLADI';  ELEMENTS 

EUIMENTAL  WAKi:  SURFACE  DESCRIPTION 


PRESCRIBED  CUSSIC  WAKE  MODEL 

NUMBER  OF  SPAiN  SEGI<IENTS  15 

NUMBER  OF  STREAI4WISE  SEGMENTS  TT 

WAri  ELEMENT  ANGULAR  IIICHEMENT 

a)  NEAR  iiLADE  REGION  20° 

b)  FAR  BLADE  REGION  (?•  90°  FROM  BLADE)  35° 

NUMBER  OF  WAKE  TURNS  COMPLETED  7 .  3 

TO'TAL  NUMBER  OF  BASIC  WAKE  ELEI-HJNTS  1155 


TABLE  3 
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COMPUTER  EXECUTION  TIME  SUMMARY 
FOR 

CAGi;  STUDY  2 


EXECUTABLE  SECTION 

I.  INFLUENCE  CALCULATIONS  PHOUHAM 

A.  BLADE  INFLUENCE  CALCULATIONS 
(TOTAL  NUMBER  =  720,000) 

ii.  WAKE  INFLUENCE  CALCULATIONS 
(TOT'AL  NUMBER  =  1,386,000) 

C .  OTHER  ROUTINES 

II.  CREATE  SOLUTION  FILE  PROGRAM 

III.  SOLVE  BY  ITERATION  PROGRAM 

(TO'TAL  NUMBER  OF  ITERATIONS  =  28l ) 

IV.  TRANSFER  SOLUTION  FILE  PROGRAM 

V.  VELOCITY  CALCULATION  PROGRAM 

TOTAL  ABOVE 


TIME 

(SECOfRJS) 

2328 

5‘»75 

20 

137 

5970 

U 

130 

ll*33l< 


TABLE  U 
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APPl'iriDIX  A 

nur^niARY  of  the  trameformation  equationh 


Preceding  page  blank 
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We  are  concerned  here  with  presenting  in  suinmai^  form  the 
transformation  equations  which  relate  an  inertial  reference  fran-.e 
(  / ,  '/  /  /  )  to  a  body  fixed  reference  frame  ^  ) 

in  which  the  body  is  translating  with  a  constant  velocity  k/ 
given  by 

Vj  -  I  (A.l) 


and  rotating  about  the  i  axis  at  a  constant  rotational  rate 
given  by 

We  assume  that  at  time  y *  /  ' '  ^  the  axis  systems  were  coincident. 
It  follows  then  that  the  coordinates  of  a  point  '  in  the  inertial 

frame  are  related  to  the  body  fixed  frame  as  follows: 


X  X  ^  ^  (A. 3) 

y'  -  /<  y  ^  y ^ic  ^  (a. It) 

-  2  -  Vj.  y  (A.  5) 

where  the  inverse  is  given  by 

y  -  X '  '  ~  y  ^  ^ '  ( A .  6 ) 

y  "  y '  ^  y  (A. 7 ) 

y  --  2'  f  \/^y'  .  (A. 8) 


Ill 


From  these  relations  it  follows  then  that 


- 

y'  ^ 

^  ^  Jy' 

(A. 9) 

2 

j^/ 

(A. 10) 

h 

2 

Ji  '  . 

(A. 11) 

Fxirthermore ,  derivations  with  respect  to  are  given  by 


/^_7'l7Vy  / 


where 


17  ’ 


'  J>'  '  J/- 


^  ^  J  J  ^ 

^  X  f  y  y  ^  ^ 


It  follows  from  the  above  that 


V  '■  ^7 


and  similarly 


V7  ^  V  ' 


(A.i:) 


(A.lM 


(A. 15) 


(A. 16) 


■ 
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The  substantial  derivations  are  related  as  follows 


-  ) 

-  '  J 

y-  V  =  jy'  ^ 

k  '  \7  '  ^ 

(A. IT) 

where 

/ 

/S 

xs 

V 

--  X 

y  y  ^ 

(A. 18) 

^  / 

t 

i7 

--  X 

/  Vy  y  ^ 

y/^  . 

(A. 19) 

The  velocities  are 

related  as 

follows 

V 

--  / 

U'  ^  /i  y  K 

/ 

^A.20) 

Furthermore,  any  scaler  function  /V  evaluated  at  some  point 

^  J  must  have  the  same  value  when  evaluated  at  point 
y  it  ‘  / 'J  if  point  is  related  to  point  ^  usin^ 

the  set  of  transformation  equations  given  by  equations  (A. 3),  (A.l*) 
and  (A.5).  That  is 

=  y//y;y;y;y'J 

(A. 21) 
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Using  relation  (A.  12)  we  may  note  that  for  some  scalar  function 

which  is  independent  of  /'  that  the  above  relation  reduces 
to 

y ! /'J  ^  (a.  22) 
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APPENDIX  B 

RELATIONS  BETWEEN  DOUBLETS  AND  VORTICITY  AND  EXTENSION 
OF  HEIWHOLTZ  CONDITION 


We  will  discuss  in  this  section  the  relations  between  doublets 


and  vorticity.  Having  established  the  relationships  it  is  then  an 
easy  matter  to  extend  Helmholtz's  theory  which  concerns  vorticity  to 
a  corollary  involving  doublets  in  order  to  determine  the  strengths 
of  the  wake  doublet  surface  distributions. 

Let  us  assume  that  we  have  a  surface  doublet  distribution  whose 
potential  at  some  point  is  given  by  /"V  and  whose  strength 

is  given  by  •  The  doublet  axis  is  assumed  to  be  directed 

everywhere  along  the  positive  local  surface  normal  direction.  We  shall 
assiune  further  that  and  are  two  orthogonal  axes  lying  on  the 

surface  with  the  X  axis  directed  along  the  unit  normal  as  shown  in 
Figure  B.l.  If  doublets  and  vorticity  are  related  we  should  be  able 
to  derive  an  expression  relating  the  local  doublet  strength  /y 
to  the  vorticity  components  and  ^ J  whose 

assumed  directions  are  as  shown  in  Figure  B.l.  In  order  to  derive  this 
relationship  we  shall  assume  that  U,  V  euid  W  are  the  velocities  along 
the  /y  and  X  axis  respectively.  It  can  be  shown  now  that 

^  (b.d 

where  X'  and  X  represent  some  small  distance  above  and  below  the 
vortex  location  XX  'X /  ' 

The  velocity  i/  ,  however,  is  related  to  the  derivatives 
of  the  doublet  potential  with  respect  to  /  such  that  we  may  write, 
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(B.2) 


If  we  now  rewrite  the  derivative  in  terms  of  its  definition  as  a  limiting 
process  we  have 


y/',  >}) 


A/ 


V 


r  'J  '  ■'),  /  J 


Af 


(B.3) 


Rearranging  the  terms  above  we  may  further  write  ^  as , 


y  ^  A  J 

'  A/ -'4. 


y),  r-'J  -  J 


Z)/ 


r'J  -  ^ 

A/ 

(B.I4) 


It  cem  be  shown  that  the  value  of  the  doublet  potential  at  a  point 
Just  above  or  below  the  doublet  itself  approaches  one  half  the  negative 
or  positive  value  of  the  Iona t  doublet  strength  itself.  Thus  we  find 


that 


AJ 


^  Ar 


'Ir^^ 


(B.5) 


But  the  limit  term  is  Just  the  definition  of  the  derivative  of  the 
doublet  strength  with  respect  to  /  ,  therefore  the  relationship 

between  the  ^  vorticity  and  doublet  strength  at  a  point  is  given  by 


(B.6) 


In  a  similar  manner  or  by  using  the  following  vortex  compatability 
relation 


b/  -  il 
-  n 


(B.T) 


it  follows  that 


(B.8) 


Now  about  a  wing  surface  it  can  be  shown  that  at  a  particular  span 
station  yy  the  total  circulation  J  is  given  by  the  line 

integral  of  the  vorticity  from  the  point  on  the  trailing  edge 

lower  surface  J  point  on  the  trailing  edge  upper 

surface  /  that  is  referring  to  Figure  B.2, 


that  is  referring  to  Figure  B.2, 
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<«//'  (B  9) 

Using  relation  B.fi  above  it  follows  then  that  the  circulation  about  a 
particular  space  station  is  given  in  terms  of  doublet  strength  by, 

^<7^^  j  -jj  J  . 

(B.IO) 

In  the  wake  region  the  application  of  Helmholtz's  theorem  on  the 
conservation  of  vorticity  shows  that 


J7 


(B.ll) 


where  now  /  and  'p  are  orthogonal  ajcej  lying  on  the  wake  surface 
with  the  aucis  directed  edong  a  local  streamline  and  the  Z 

axis  being  directed  generally  from  root  to  tip  across  the  wake.  In 
order  for  a  zero  pressure  discontinuity  to  exist  across  the  wake  in 
steady  flows  it  can  be  shown  that  nmst  be  zero.  Thus  it  follows 

using  relation  (B.6),  where  we  now  define  as  the  local  wake 

surface  doi  t  strength  whose  axis  is  along  the  Z  axis ,  that 


(B.12) 


everywhere  in  the  wake  surface.  Since  Kutta's  hypothesis  requires 
the  fluid  to  leave  the  trailing  edge  smoothly  and  because  the  above 
relation  implies  that  is  a  constant  along  a  streamline  we 

find  following  such  a  streamline  that 

.  (B.13) 

Thus  we  find  that  along  a  wake  streeunline  the  value  of  the  weke  doublet 
strength  is  a  constant  and  is  given  by  its  value  at  the  trailing  edge, 
that  is 

along  a  wake  streamline. 

(B.iM 

Note  at  this  point  that  by  properly  integrating  about  the  wake  region 
one  finds  (as  one  should)  that  the  circulation  in  the  wake 
about  some  constant  to  be 

which  is  Just  the  negative  of  the  bound  circulation  given  by  equation 
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APPENDIX  C 


ANALYTIC  EXPRESSIONS  FOR  AND 
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We  vrlll  be  concerned  in  this  section  with  deriving  an  analytic 
algebraic  expression  for  the  doublet  potential  and  its  derivative 

^  ^  at  some  field  point  given  some  eurbitrary  >i  sided 

planar  doublet  distribution  of  constant  strength  whose  doublet  axis  is 
everywhere  normal  to  the  surface.  Before  we  begin  this  derivation 
let  us  first  define  the  doublet  so  as  to  make  clear  the  axis  convention 
used. 

In  a  physical  sense  a  doublet  is  derived  as  the  sum  of  a  source 
and  a  sink.  Suppose  we  have  a  point  source  of  strength  ^  and 
a  point  sink  of  strength  located  about  the  origin  of  an  y,  y  J 

rectanguleur  coordinate  system  as  shown  in  Figure  C.l.  It  follows  from 
this  that  the  potential  at  some  field  point  P  is  given  by 


(C.l) 


where  the  position  vectors  and  /<'.?  are  as  shown  in  the  figure. 

Let  us  now  define  the  vector  distance  from  the  sink  to  the  source  as 
^  where 


(C.2) 


where  k  is  the  unit  vector  along  the  line  Joining  the  sink  to  the 
source.  It  follows  from  this  that  the  doublet  potential  at  the  field 
point  P  is  obtained  by  taking  the  limit  of  the  source  plus  sink 
potentieil  at  the  point  P  as  they  each  approach  the  origin  assuming 
the  product  of  the  strength  and  vector  distance  ^  remain  constant 


AXIIJ  CONVENTION  RELATING  TliE 
POTENTIAL  AT  A  POINT  INDUCED 
iiY  A  GOUHCE  +  GINK 


FIGURE  C.l 
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auid  given  by  where 


--  y*  -  ^  M  ,  (C.3) 

Performing  the  limiting  process  we  find  that  the  potential  at  a 
field  point  resulting  from  a  doublet  situated  at  the  origin 

in  Figure  C.l  is  given  by 
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where  is  the  direction  of  the  doublet  axis  directed  along  the 

positive  /  axis.  Extending  this  result  to  a  surface  doublet 
distribution  whose  axis  is  everywhere  directed  normal  to  the  local  sur¬ 
face  we  find  that  the  potential  at  a  point  becomes 

/  r  My/  • 

'  '  ^7y  /  ■  ■■  ;  -  /  y~ 

where  /(j,--  is  the  vector  distance  from  a  point  on  the  local  doublet 
surface  to  the  field  point  and  is  now  the  doublet  strength 

per  unit  area. 

Given  this  expression  for  the  potentieLl  what  we  now  wish  to  do  is 
to  integrate  this  expression  for  an  arbitrary  >7  sided  plemar  doublet 
surface  distribution  of  constant  strength.  More  specifically  we  may 
state  the  problem  as  follows: 


Given:  A  planar  ^  sided  figure  as  shown  in  Figure  C.2  (for 
H  =  U)  whose  surface  lies  within  the 
plane  and  whose  positive  unit  normal  ?t  is  along  th^ 
positive  y  axis  and,  given  this  plane  to  be  a  surface 
doublet  distribution  plane  of  constant  strength 
whose  axis  is  everywhere  directed  along  the  unit  normal 
then 

Determine:  A.  The  potential  at  an  arbitral^  point  ^ 

whose  coordinates  are  /X,  y,  ^  J 
B.  The  vector  gradient  of  ^  at  some  arbitrary 
point  .that  is, determine 

where 


(c.f.) 
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Note,  in  the  analysis  to  follow  the  point  coordinates  are 

given  in  terms  of  an  reference  system  which  is  coincident 

with  the  y J  system  (refer  to  Figure  C.2)  in  order  to  distinguish 

the  coordinates  of  the  fixed  point  from  the  surface  integration 

variables. 

Now  from  Figure  C.2  we  may  identify  the  following  relationships: 

^  r  (c.T) 

''  ^  A  ^  (C.8) 


AXII3  CONVENTION  RELATING  Tlll'i 
PLANAR  N=U  GIUEU  DOUBLET  SURKACE 
DISTRIBUTION  TO  GOME  FIELD  POIiJT 


FIGURE  C.2 
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The  problem  now  is  the  evaluation  of  the  integral  ^ .  Given 
the  evaluation  of  Z' ^  then  the  gradient  may  be  obtained  by 

straightforward  differentiation.  Note  that  the  gradient  may  be  taken 
inside  the  integral  and  then  the  resulting  integrand  may  be  integrated. 

We  choose  not  to  use  this  procedure  as  singtdarities  of  higher  order 
than  what  already  exist  will  result.  Althougli  these  singularities  may 
be  evaluated  in  terms  of  the  Cauchy  principal  values  it  does  unnecesseu’ily 
complicate  the  evaluation. 

Let  us  consider  now  the  evaluation  of  the  integral  given  by 


(C.13) 


In  Figiu*e  C.3a  we  have  sketched  the  surface  S  for  an  X  =  L 
sided  figure  and  have  indexed  the  corner  points  in  a  manner  such  that 
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INDEXING  CONVENTION  FOR  THE 
PLANAR  N=l»  GIDED  GURFACE 


FIGURE  C.3a 

GUbDlVIGlON  OF  REGION  S  INTO  HEGIONG 


FlGUl^E  C.3b 


when  viewing  the  surface  adong  the  negative  surface  normal  direction 
the  corner  points  are  indexed  consecutively  ^  in  a 

clockwise  direction.  Figure  C.3b  shows  that  the  surface  S  may 
be  subdivided  into  H  regions  such  that  the  integral  may  be  evaluated 
as  the  sum  of  Integrals  each  integrated  about  a  specific  region, 

.  We  may  define  these  regions  more  precisely  if  we  first 
let  the  corner  point  coordinates  be  given  by  ^ define 
J"  as  the  indexed  corner  point  located  immediately  clockwise  (with 
respect  to  the  negative  surface  normal  direction)  to  the  indexed 

corner  point.  Now  the  equation  of  the  straight  line  Joining  corners  / 
and  J  is  given  by  j  -  where 

^  ^  ^  (C.U) 

where 

/  lies  between  /  and  /' 

^  -K  J  . 

Note  at  this  point  that  if  the  slope  V  straight  line 

is  infinite  a  potential  problem  exists.  However,  in  this  case 

-  /-  and  the  surface  aurea  of  the  region  j  is  zero.  Thus 

the  integral  contribution  of  this  region  is  also  zero. 
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It  follows  now  that  the  integral  may  be  written  as 


^  /  y  //^.'y)  (c.iT) 


where 
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V  y  ^  V 
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(C.18) 


The  integration  with  respect  to  may  however  be  rewritten  as 

follows 


* 
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(C.19) 


Thus  we  may  write 


(C.20) 


The  second  summation  term  is  zero,  however,  since  after  integrating 
with  respect  to  ^  the  integration  with  respect  to  Z'  effectively 
represents  a  line  integral  evaluation  about  a  closed  path.  Thus  the 
integration  reduces  to 


Z' 


(c.?i) 


where  we  have  defined 


(C.22) 


The  integration  with  respect  to  may  be  carried  out  using 

Reference  15  such  that  the  resulting  expression  for  becomes 


- 

(C.23) 

The  evaluation  of  this  remaining  integral  may  be  carried  out  using 
Reference  16  after  substituting  for  ^  ^  using  relation  (C.II4). 
The  resulting  analytical  algebraic  relation  for  ^ Z/^J  is  given  by 
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Summarizing  what  we  have  done  to  this  point,  we  can  write  the 
analytic  algebraic  expression  for  the  doublet  surface  potential  at 
some  field  point  as 


where 
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(C.P5) 


(C.26) 


(C.2T) 

(C.28) 


(C.20) 

(C.30) 

(C.31) 


(C.32) 


The  cartesian  gradient  vector  of  the  potential  can  now  be 


determined  by  simply  differentiating  equation  (C.P5).  The  result  may 
be  wi-itten  as 


y 


/ 


(c.r' 


where  for  5  a  dummy  variable  of  differentiation  we  find 


_ 
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where 


/  ^  r  ' 


(C.3^) 
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I  J3  ^  '■  (C.37) 


Z/  /.  --  ^  (C.38) 

/  \J  /  * 


The  renaining  parameters  are  given  by  equation  (C.28)  through  (C.3l). 

The  expression  for  is  valid  for  all  Pj 

locations  except  for  a  point  ^  lying  on  the  edge  of  the  surface 
where  it  is  undefined.  For  a  point  in  the  plane  of  the  surface 

we  find  that  and  are  zero. 

The  expression  for  is  valid  for  all  >;  yJ 

locations  away  from  the  doublet  plane.  For  the  point  lying  in 

the  plane  of  the  doublet  distribution  we  find  that 

A.  For  within  the  doublet  surface 

y  *  (C.30) 

B.  For  outside  the  doublet  surface 

~  (C.I40) 

C .  For  X  on  the  edge  of  the  doublet  surface 

/X,  VyPj  ~  UNDEFINED.  (C.l*l) 

We  should  note  at  this  point  that  in  the  computer  program  all 
velocities  and  lengths  are  nondimensionalized  based  on  the  rotor  tip 
speed  X'/x>X)  and  the  rotor  radius  respectively.  Hince 

actually  represents  the  velocity  Xs  along  the  S 
direction  at  point  Z’  the  form  of  the  preceding  equations  as 
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propranmed  is  aa  follows  where  the  tllda  ^ J  represents  the 

iiondimensional  parameter: 

.  -f  <' 


where 
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(C.U3) 
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(C.I4U) 

^  " 

(C.I15) 

The  actual  < 

expression  for  ^  is  the  same  as 

for 

A 

but  now 

all  lengths 

are  to  be 

first  nondimensionalized  by 

the 

radius 

A  , 

for  example 

instead  of 

/  read  ?  where 

'Xd 

R 

This  completes  the  discussion  of  the  analytic  expressions  for 
and  //>J  . 
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APPENDIX  D 

DETAILED  DISCUSSION  OF  THE  COMPUTER  PROGRAM 
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D.l  INTRODUCTION 

In  the  following  discussions  we  will  consider  the  internal 
structure  of  the  program  at  the  level  of  the  subroutines  in  termr.  of 
their  function.  Following  this  we  will  discuss  the  program  operation 
and  where  necesseiry  a  specific  flow  diagram  will  be  presented.  The 
above  discussion  ajid  flow  diagram  will  then  serve  as  an  aid  in  describing 
the  major  program  options  and  the  input  and  output  procedures. 

In  Section  U.O  we  noted  that  the  program  is  divided  into  seven 
major  programs  each  of  which  is  concerned  with  a  particular  aspect  of 
the  overall  problem.  In  the  discussions  to  follow  we  will  consider 
these  programs  one  at  a  time  and  discuss  the  programs  as  though  they 
were  independent  of  each  other.  It  will  be  assumed  here  that  the 
reader  is  familiar  with  the  symbols  and  coordinate  reference  syntemn 
described  earlier  in  Sections  (?)  and  (3). 
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D,?  t^AIN-INFLUENCE  CALClILATIONn  PROORAM 

In  this  section  we  will  discuss  the  program  MAIN  -  INFLUENEK 
CALCULATIONS.  It  is  the  function  of  this  program  to  define  the 
necessary  geometry  and  to  perfoim  the  calculations  necessary  to  define 
the  blade  and  wake  influence  coefficients.  We  shall  present  first  a 
summary  discussion  of  each  subroutine  found  in  this  program. 
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D.2.1  Suraraary  Discussion  of  the  Subroutines 


The  subroutines  used  in  the  MAIN  -  INFLIJENCK  f'ALCUIJVTIorK’  prorram 
are  summarized  here  in  the  general  order  in  which  they  are  railed 
upon  b-,  the  main  program. 

Subroutine  irEAD(K) 

The  function  of  this  subroutine  is  to  simply  print  a  heading  page 
indicating  the  main  program  title.  This  page  is  used  as  a  cover  page 
for  the  results  printed  out  later. 

Subroutine  TNI 

The  function  of  this  subroutine  is  to  accept  the  card-itiput 
necessary  to  define  the  rotor  blade  for  the  specialized  planform 
option  described  in  section  3.2.  The  routine  also  checks  certain  input 
before  allowing  the  program  to  continue  to  the  next  step.  If  an  input 
error  is  found  the  program  will  abort  after  printing  out  the  data  as 
inputted  and  additional  parameters  determined  on  the  basis  of  innut 
supplied  which  would  have  been  passed  on  to  other  routines. 

Subroutine  IN2 

The  function  of  this  subroutine  is  to  accept  the  card  input 
necessary  to  define  the  rotor  blade  for  the  general  planform  option 
described  in  section  3.2.  It  also  defines  certain  parameters  based  on 
the  input  supplied  and  passes  these  on  to  other  routines. 


. . .  ■  ■■- 
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Subroutine  IN 3 

The  function  of  this  subroutine  is  to  accept  the  card  input 
necessary  to  define  the  rotor  free  stream  conditions ,  the  number  of 
rotor  blades  and  parameters  necessary  to  define  the  method  in  which 
the  wake  influence  calculations  are  to  be  performed. 

Subroutine  INH 

This  is  an  auxiliary  subroutine  which  is  not  used  in  the  rotor 
problem.  It  is  an  input  subroutine  compatable  with  IN3  and  is  used 
to  input  card  data  comparable  to  IN3  when  the  program  is  selected  to 
do  a  lifting  surface  theory  analysis  of  a  planar  wing  rather  than  of  a 
rotor  system. 

S'Ubroutine  OlJTl 

The  function  of  this  subroutine  is  to  print  the  input  data  of 
nil.  This  output  serves  as  a  permanent  record  of  the  specialized 
planform  geometry  specified  for  the  rotor  system. 

Subroutine  OUTP 

The  function  of  this  subroutine  is  to  print  the  input  data  of 
iriP.  This  output  serves  as  a  permanent  record  of  the  general  nlanform 
geometry  specified  for  the  rotor  system. 

ljubroutine  OUT  3 

The  function  of  this  subroutine  is  to  print  the  input  data  of  I'I3. 
This  output  serves  as  a  permanent  record  of  the  specified  free  stream 
condition  and  wake  calculation  mode  selected. 
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Subroutine  OUTU 


The  function  of  this  subroutine  is  to  print  the  input  data  of 
INU.  The  output  is  comparable  to  0UT3  except  that  the  output  is  for 
the  planar  wing  rather  than  the  rotor  system. 

Subroutine  OUT A 

The  function  of  this  subroutine  is  to  print  out  a  coordinate 
description  of  the  airfoil  chosen  in  the  specialized  planform  option. 
This  may  be  used  to  check  the  accuracy  with  which  the  airfoil  is 
described  in  the  program.  This  subroutine  supports  ROTORG  and  calls 
on  AFOIL. 


Subroutine  AFOIL  (X ,7.U,ZL,DZUDX ,DZLDX,BU,BL) 

The  function  of  this  subroutine  meiy  be  stated  as  follow:;; 

Given:  A  chord  station  (X),  where  X  represents  the  chord  non- 

dlmenslonallzed  distance  from  the  airfoil  leading  edge  along 
the  chord  line. 

Find:  At  the  station  X,  the  nondimensionalized  upper  and  lower  airfoil 

thickness  (ZU  and  ZL),  the  upper  and  lower  airfoil  tangent 
slopes  (DZUDX  and  DZLDX)  and  the  angle  the  airfoil  tangent  lines 
make  with  the  chord  line  in  radians  (BU  and  BL). 

The  subroutine  defines  the  NACA  OOXX  family  of  airfoils  as  a  set  of 
equations.  The  thickness  ratio  (TC),a  parameter  in  this  equation, is 
transferred  to  the  subroutine  implicitly.  This  subroutine  supports 
ROTORG , 
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Subroutine  ARCL  (XI ,XF,DARCLU ,DARCLL) 

The  function  of  this  subroutine  is  to  cedculate  the  upper  and  lower 
airfoil  arc  distances  (DARCLU  emd  DARCLL)  lying  between  two  lines  drawn 
perpendicular  to  the  two  chord  stations  (XI, XF).  This  subroutine  calls 
on  AFOIL  and  is  generally  not  used  in  the  overall  problem.  Tt  serves  a 
part  in  checking  the  planform  description  when  this  option  to  be 
described  later  is  chosen.  This  subroutine  supports  ROTORC. 

Subroutine  GEOMRl  (XC,YTS,ITIP) 

The  function  of  this  subroutine  may  be  stated  as  follows: 

Given;  A  control  point  location  in  terms  of  a  chord  station  (XC) 
and  a  true  span  station  (YTS). 

Find:  The  actual  upper  and  lower  surface  coordinates  and  the  local 

surface  tajigent  slopes  in  terms  of  the  blade  one  coordinate 
system. 

The  information  is  transferred  implicitly.  The  parameter  ITIP  is  used 
to  indicate  whether  the  rotor  tip  or  tips  are  to  be  faired.  This 
subroutine  supports  ROTORG. 

Subroutine  RO'l'ORG 

The  function  of  this  subroutine  is  to  serve  as  a  control  nrorram 
in  order  to  systematically  define  the  actual  control  point  locations  in 
terms  of  the  blade  fixed  coordinate  system  when  the  specialized  planform 
description  option  is  selected.  Its  main  supporting  subroutine  is 
OKOI'iRl. 


Ill  aitdition  to  this  function  the  subroutine  has  an  option  which 


nerves  the  purpose  of  printing  the  surface  grid  point  coordinates 
and/or  card  punching  these  coordinates  in  a  form  compatable  with  that 
required  as  input  data  by  subroutine  IN2. 

I'lirthermore ,  this  subroutine  contains  an  option  which  will  calcu¬ 
late  mean  control  surface  locations  and  slopes  and  print  these  results. 
'I'his  option  may,  however,  not  be  selected  during  the  running  of  the 
overall  program.  Independently  this  option  serves  as  a  reference  check 
on  the  elemental  control  surface  calculations  done  in  a  different  manner 
later  in  the  program. 

Subroutine  SURF 

The  function  of  this  subroutine  is  to  define  an  elemental  planar 
four  sided  control  surface  in  terms  of  its  location,  corner  point 
coordinates  and  transformation  matrix  given  initially  four  control 
point  locations  which  do  not  necessarily  lie  in  a  plane.  This  sub¬ 
routine  supports  CCURF  and  CWAKE. 

Subroutine  cn  (CX,Y,XB,YB,ASUM ,  A) 

The  function  of  this  subroutine  is  to  calculate  the  centroid 
location  (X,Y)g,  and  area  (ASUI'l)  of  a  four  sided  planar  figure  whose 
corner  point  coordinates  are  given  by  (X,Y)j^_2  ^  where  (X,Y)^  =  (0,0) 
and  the  corners  are  numbered  in  a  clockwise  direction.  This  subroutine 
supports  OURF ,  and  APPROX  . 
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r;ubroutine  ROTZ  (Tir.TA,X,Y,AJ) 

The  function  of  this  subroutine  may  be  stated  as  follows: 

Given:  A  coordinate  system  //  ,  T )  whose  origin  and  transformation 

matrix  eure  given  with  respect  to  an  {Y^Z)  coordinate  system  by 
X,Y,Z  and  AJ. 

Find:  The  origin  and  transformation  matrix  of  the  system 

if  it  is  simply  rotated  about  the  Z  eocis  through  a  displacement 
THETA.  THETA  being  positive  in  the  usual  sense. 

The  results  are  transferred  explicitly  using  X,Y  and  AJ.  This  subroutine 
supports  CWAKE  and  CBLADE. 

Subroutine  CALLOC 

The  function  of  this  subroutine  is  to  reallocate  the  main  core 
storage  registers. 

.Subroutine  REMUM 

The  function  of  this  subroutine  is  to  renumber  the  blade  elemental 
control  surfaces  in  terms  of  an  integer  I.  Prior  to  the  use  of  this 
subroutine  the  upper  ajid  lower  control  surfaces  were  sequenced  in  terms 
of  (N,M). 

Subroutine  WTl  (ITI.IRW) 

The  function  of  this  subroutine  is  to  write  or  read  off  external 
file  ITl  (depending  on  whether  IRW  equals  0  or  is  not  equal  to  0 
respectively)  the  control  parameters  and  data  allocated  in  designated 
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common  storage  locations  LINDA  and  BRENDA.  Dome  parameters  and  data 
in  these  two  locations  are  needed  effectively  in  all  of  the  major 
programs.  Thus  this  subroutine  and  file  ITl  serves  as  a  continuity 
link  among  the  programs. 

Subroutine  CSURF 

The  function  of  this  subroutine  is  to  serve  as  a  control  program 
in  order  to  systematically  define  the  geometry  of  all  the  blade  one 
elemental  control  surfaces.  The  main  supporting  subroutine  used  for 
this  task  is  SURF. 

In  addition  this  subroutine  computes  the  free  stream  veloci  ty 
conditions  on  each  blade  element  and  writes  all  the  above  results  on 
external  file  ITl. 

Subroutine  INTEG  (X,Y,Z,G,E,XI,XIS) 

The  function  of  this  subroutine  is  to  calculate  the  nondimen- 
sionalized  doublet  potential  and  velocity  influence  coefficients  in 
/J  reference  system  at  a  field  point  ,Y,Z)  resulting 
from  an  M  sided  planar  doublet  distribution  whose  direction  is  along 
the  T  axis  and  whose  corner  point  locations  are  given  clockwise 
in  the  //  '^j  plane  as  (G,E)[^^.  The  potential  influence  coefficient 
is  given  by  XI,  and  the  velocity  influence  coefficients  in  the  /,  y 
and  f  directions  eire  given  by  XIS(l),  XIS{2)  and  XIS(3)  respectively. 
The  subroutine  tests  to  ensure  that  the  field  point  does  not  lie  on 
the  edge  of  the  planar  surface.  This  subroutine  supports  CBLADE  and 


CWAKE  subroutines. 
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LJubroutlne  CBLADl']  ( IGTAHT.IGTOP) 

The  function  of  this  subroutine  is  to  serve  as  a  control  profiram 
in  order  to  systematically  define  the  total  influence  coefficients 
of  all  J=l,Ji<lAX  blade  elements  on  blade  one  elements  I  where  I 
lies  inclusively  between  control  input  specified  elements  l=lGTAiH' 
and  I=IGT0P. 

The  blade  influence  coefficients  are  calculated  with  respect  to 
the  (X' ,Y' ,Z' )  blade  fixed  axis  system  and  are  written  out  on  file 
IT2  immediately  after  all  J  elements  of  all  blades  influence  on  a 
blade  one  element  I  are  kncwr,  '^hey  are  written  out  as  three  sets  of 
JMAX  influences  corresponding  to  tht  X',Y'  and  Z'  directional  influences. 

Subroutine  PWAKT 

'.'he  function  of  this  subroutine  is  to  serve  as  n  preliminary 
subroutine  to  the  wake  subroutine  when  specific  program  options  are 
selected.  Basically  it  is  used  when  the  blade  influence  coefficients 
were  defined  during  a  previous  execution  run  and  entry  is  now  made 
into  the  CWAKE  subroutine  directly.  It  serves  the  function  of  readin/' 
external  file  ITl  and  positioning  in  main  core  storage  the  data 
necessary  to  execute  the  CWAKE  routine. 

Subroutine  CWAKE 

The  function  of  this  subroutine  is  to  serve  as  a  control  program 
in  order  to  systematically  define  the  total  wake  influence  coefficients 
of  all  M=1,MM  wake  span  stations  on  blade  one  elements  1=1,  PlAX . 


The  procediire  is  to  begin  at  span  station  M=1  and  to  proceed 
streamline  wise  down  the  wakes  calculating  veike  elements  and  their 
influence  coefficients  with  respect  to  the  {X',y',7.')  system  on  all 
1=1,  IMAX  blade  elements.  At  any  single  wedte  span  station  the 
influences  are  additive  thus  a  single  set  of  IMAX  total  influence 
coefficients  is  kept.  When  the  final  streamwise  waJie  element  trailinj' 
M=1  wake  span  station  is  defined  the  set  of  three  vectored  IMAX  total 
wake  Influence  coefficients  are  written  out  cn  external  file  IT3.  The 
procedure  is  now  repeated  for  span  station  M=M+1  through 

This  subroutine  also  has  a  planar  wake  model  to  be  used  during 
the  planar  wing  option. 

oubroutlne  SETUPl 

The  function  of  this  subroutine  is  to  Interchange  the  IMAX 
column  and  three  sets  of  MM  rows  of  data  written  on  external  file 
IT3  during  the  CWAICE  routine  and  write  this  information  on  file  ITl 
behind  the  information  already  written  on  ITl. 

Subroutine  A^^PROX  ( X ,Y  ,Z ,r.  ,E,XT  ,XTG ) 

The  function  of  this  subroutine  is  to  replace  TNTES  to  calculate 
the  influence  coefficients  when  the  field  point  (X,Y,Z)  is  farther 
than  a  certain  distance  from  the  centroid  of  the  doublet  element. 

It  uses  multi-pole  expansion,  which  is  an  approximation  to  the 
algorithm  shown  in  Appendix  C,  and  will  save  the  time  of  calculation. 
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D.2.2  Program  Operation  in  Terms  of  a  Flow  Diagram 

Having  completed  now  the  discussion  of  the  subroutines  we  may 
now  describe  the  "Influence  Calculation"  program  in  terms  of  a  flow 
diagram.  In  Fig\ire  D.l  is  shown  the  detailed  flow  diagram.  We  have 
attempted  to  present  this  flow  diagram  at  the  level  of  the  sub¬ 
routines.  We  tiave  indicated  in  the  flow  diagram  the  stage'"  at  which 
input  is  required  and  where  output  is  effected.  We  have  shown  on  the 
side  of  the  flow  diagram  proper  the  supporting  elements  including 
external  files  and  supporting  subroutines  required  by  each  major  sub¬ 
routine  called  on  by  the  major  program. 

Note  in  this  diagram  that  the  program  has  essentially  5  entry 
points  and  5  exit  points  depending  on  the  value  of  the  parameters 
MCTNL(l)  and  MCTRL(P)  respectively.  Note  further  that  the  parameter 
MCTIU-(')  determines  whether  INI  (Rotor  Analysis)  or  Illlj  (Planar  Wing 
Analysis)  is  selected.  The  external  files  are  designated  as  I'i'l, 

ITP  and  IT3.  If  the  program  is  executed  from  beginning  to  end  in  one 
step  three  files  are  required.  Note,  however,  that  by  executing  the 
program  in  stages  using  the  program  option  parameters  MCTRL(l)  and 
MCTRLCp)  the  maximum  number  of  on  line  external  files  needed  at  it.y 
one  time  is  two  if  the  nur.iber  of  oiade  elements  is  1000  or  less,  i'he 
third  file  in  this  case  must  be  defined  as  a  dummy  file,  say  equal  to 
one  of  the  other  two.  Files  ITl  and  IT2  on  program  completion  contain 
all  the  information  needed  in  later  main  programs  and  sire  to  be 
considered  as  permanent  storage  files.  File  IT3  on  program  completion 
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nil  filtered  format  on  file  ITl.  As  such  Kile  IT3  may  be  dlscar'leci  on 
c-ompletion  of  the  MAIN  -  INFUJKNCK  CALCULATION;-.  PROGRAM. 

1).?.  R  Program  Input  Procedures 

The  input  required  and  input  sequence  can  be  determined  from  the 
flow  diaRram  and  this  of  course  is  a  function  of  the  MCTRL  parameters 
which  are  inputted  in  the  NAMELIST/MAINl/  data  group  at  the  start 
of  the  program. 

Lxcept  for  the  input  read  in  subroutine  IN?  and  the  identification 
name  assigned  to  the  run,  the  input  is  entered  in  a  NAMKI.IUT  mode. 

This  allows  us  to  preset  the  input  parameters  and  only  input  those 
parameters  which  we  want  to  override  in  any  particular  run. 

We  shall  consider  now  the  details  of  inputting  the  data  in  the 
order  in  which  they  ore  executed  in  the  overall  program,  that  is 

1.  Main  Program  Control  Parameters 

NAMELIGT/MAINl/ 

2.  Message  Card 

NAME 

3.  a  Specialized  Planform  Geometry  Description 

m ELIGT/DATA?/ 

or 

3.b  General  Planform  Geometry  Description 

ti.a  Rotor  Free  Stream  and  Wake  Analysis  Parameters 
NAI-IELIST/DATA3/ 

or 

li.b  Planar  Free  Stream  and  Wake  Analysis  Parameters 

NAMELIST/UATAU/ 


ir.o 

1.  Main  Program  {,'ontrol  Parameters 

The  following  data  is  to  be  inputted  in  the  NAMELinT/fiAIWl/ 

data  statement: 

MCTHh(l)  =  K  (Default  value  is  K  =  1) 

K  =  5  if  the  program  is  to  effect  execution  at  GETUPl. 

K  =  !<  if  the  program  is  to  effect  execution  at  CWAKE. 

K  =3  if  the  program  is  to  effect  execution  at  CBLADE. 

K  =  2  if  the  program  is  to  effect  execution  at  the  initial 

entry  point  and  use  the  General  Planform  Option. 

K  =  1  if  the  program  is  to  effect  execution  at  the  initial 
entry  point  and  use  the  Specialized  Planform  Option 

MC‘rdL(2)  =  K  (Default  value  is  K  =  0 

K  =  It  if  the  program  is  to  terminate  execution  after  compietion 
of  CWAKE. 

K  =  3  If  the  program  is  to  terminate  execution  after  completion 
of  CblADE. 

K  =  2  if  the  program  is  to  terminate  execution  after  completion 
of  CSURF. 

K  =  1  if  the  program  is  to  terminate  execution  after  comf)letion 
of  HOROKG. 

K  =  0  if  the  program  is  to  terminate  execution  after  completion 
of  the  entire  program. 

MCTRL(3)  =  K  (Default  value  is  K  =  l) 

K  =  2  if  the  program  is  to  execute  the  "Planar  Wing  Option". 

K  =  1  if  the  program  is  to  execute  the  "Kotor  Blade  Option  ’ 

MCTAL(ll)  =  K  (Default  value  is  K  =  1) 

K  =  I  where  I  is  an  integer  indicating  the  blade  element  I  on 
which  the  blade  influence  calculations  of  CBLADE  are  to 
begin,  (l  I  ^  Iraax) 
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MCTRL(12)  =  K  (Default  value  is  K  =  Imax) 

K  »  I  where  I  ia  an  Integer  indicating  the  blade  element  I 
on  which  the  blade  influence  calculations  of  CBLADE 
are  to  end.  (1  '  I  ^max^ 

ITl  =  K  (Default  value  is  K  =  l) 

K  =  I  wh^re  I  is  aji  integer  indicating  the  external  file 
asa..gned  number  on  which  the  control  parameters, 
selected  geometry,  blade  element  geometry  and  on 
program  completion  wake  influence  coefficients  ai'e  to 
be  written. 


IT2  =  K  (Default  value  is  K  =  2) 

K  =  I  where  I  is  an  integer  indicating  the  external  file  on 
which  the  Blade  Influence  Coefficients  eu*e  to  be 
written. 


IT3  =  a  (Default  value  is  K  =  3) 

K  =  I  where  I  is  an  integer  indicating  the  external  file  on 
which  the  Wake  Influence  Coefficients  are  initially 
to  be  written  prior  to  entering  GETUPl. 

DUMMY  =  K 

K  =  I  where  I  is  an  integer.  If  the  default  values  of  the 
above  parameters  are  to  be  selected  simply  input  this 
parameter. 


2.  Message  Card 

Any  message  or  title  which  is  to  be  assigned  to  the  run  for 
identification  is  to  be  inputted  here.  The  message  must  be  on  one 
card  and  may  fill  part  or  all  of  the  80  column  c£U*d. 


3. a  Cpecialized  Planforra  Description 


The  following  parameters  must  be  Inputted  as  part  of  the 
IJAMELIST/DATAx/  data  statements; 

CGU  =  0.0,  A. a .  B.b,1.0  (Default  values  are  O.O) 

where  SG  is  maximally  dimensioned  for  52  real  numbers  each 
indicating  the  true  span  station  along  which  the  upper  and 
lower  surface  grid  control  points  are  to  be  defined.  The  first 
and  last  station  must  be  0.0  and  1.0  respectively,  where  f'G  = 

0.0  is  the  root  span  station  and  GG  =  1.0  is  the  tip  span  station. 
The  list  must  be  of  increasing  magnitude. 


CGU  =  0.0,  \.a  .  B.b,1.0  (Default  values  are  0.0) 

where  CGU  is  maximally  dimensioned  for  52  real  numbers  each 
indicating  the  chord  station  along  which  the  upper  surface 
grid  control  points  are  to  be  defined.  The  first  and  last  station 
must  be  0.0  and  1.0  respectively,  where  CGU  =  0.0  is  at  the 
leading  edge  and  CGU  =  1.0  is  at  the  trailing  edge.  The  list 
must  be  of  increasing  magnitude. 


CGL  =  0.0,  A.  a  .  U.h,  1.0  (Default  values  are  CGU ) 

where  CGL  is  the  same  as  CGU  except  CGL  refers  to  lower  surface 
grid  control  points.  If  no  CGL  stations  eire  Included  in  the 
/DATAl/  list  the  program  assumes  the  chord  station  grid  points 
are  to  be  symmetric  with  respect  to  upper  and  lower  surface  anc. 
thus  sets  CGL  =  CGU  for  all  CGU  input  points. 


TC  =  O.a  (Default  value  is  O.O) 

where  TC  is  the  maximum  airfoil  thickness  to  chord  ratio  of  the 
NACA  OOXX  family  of  airfoils. 


C  =  A. a  (Default  value  is  O.O) 

where  C  is  the  chord  length  to  rotor  radius  (R)  ratio. 


B1  =  O.a  (Default  value  is  0.0) 

where  B1  is  the  distance  the  root  chord  is  displaced  from  the 
axis  of  rotation  as  measiu-ed  along  the  span  axis  and  nondimen- 
sionalized  by  the  rotor  radius.  This  is  the  hub  radius  and  must 
be  greater  than  zero.  The  true  speui  is  defined  from  this  parameter 
as  (1.0  -  Bl). 


B3  =  O.a 


(Default  value  is  O.O) 
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where  133  is  the  distance  the  span  axis  is  displaced  from  the 
leading  edge  as  measured  along  the  chord  line  and  non- 
liimensionalized  by  the  chord  length. 

TO  =  A. a  (Default  value  is  O.i 

where  TO  is  the  root  airfoil  section  geometric  ang/  of  ul  taci; 
in  (leg.reeu.  '.'he  angle  is  jjositive  in  the  usual  ae"  ^dynomii* 
sense  and  is  the  angle  measured  from  the  plane  of  rotation  to 
the  chord  line. 

TT  =  A. a  (Default  value  is  O.O) 

where  TT  is  the  amount  of  lineeu*  twist  in  degrees  which  is  to 
be  applied  to  the  rotor  blade  along  the  span  axis.  It  is 
defined  as  the  tip  chord  geometric  pitch  less  the  root  chord 
geometric  pitch. 

ITIP  =  K  (Default  value  is  K  =  0) 

where  ITIO  1c  an  option  which  allows  for  fairing  of  the  rotor 
blade  ends.  Select  K  as  follows: 

K  =  P  INBOARD  TIP  AND  OUTBOARD  TIP  ID  FAIRED 
K  =  1  OUTBOARD  TIP  IP  FAIRED 
li  =  0  NEITID'JR  TIP  IP  FAIRED 


IPRINT  =  u  (Default  value  is  K  =  0) 

where  IPRINT  is  an  option  which  allows  for  printing  of  the 
computed  surface  grid  point  coordinates.  Delect  K  as  follows: 

K  =  1  Coordinates  are  printed 
K  =  0  Coordinates  are  not  printed 

IPUNCH  =  K  (Default  value  is  =  0) 

where  IPUNCIi  is  an  option  'rhich  allows  for  card  punching  of  the 
computed  surface  grid  point  coordinates  in  a  form  compatable 
for  input  in  the  "General  Planforra  Option".  Select  K  as  follows: 

K  =  1  Coordinates  are  punched 
K  =  0  Coordinates  are  not  punched 


IFOILD  »  K 


(Default  value  la  K  =  O) 


if.J) 


where  IFOILD  Is  an  option  which  allows  for  computinr  and  printing 
of  the  airfoil  coordinates  at  every  l/lOOth  chord  station. 

Delect  K  as  follows: 

K  =  1  Coordinates  printed 
K  =  0  Coordinates  not  printed 


IDURF  -  ii  (Default  value  is  K  =  0) 

where  I.^URF  is  an  option  which  allows  for  additional  elemental 
surface  parameters  to  be  computed  and  printed.  For  each  element 
an  avereige  location,  average  slope  and  calculated  surface  area 
are  presented.  If  this  option  is  selected  the  program  will 
automatically  terminate  execution  after  completion  of  this  step. 
Delect  K  as  follows: 

K  =  1  Additional  surface  data  presented 
K  =  0  No  additional  surface  data  presented 


Please  note  in  addition  to  the  above  that  the  maximum  number  of 
elemental  surfaces  that  may  be  defined  may  not  exceed  POOO.  Thus 
if  [JMAXU  and  NMAXL,  are  the  number  of  specified  chord  grid  stations 
on  the  upper  and  lower  surface  and  MHAX  is  the  number  of  sp'eci  fied 
span  grid  stations  then  the  following  relation  must  hold 

(NMAXU-1)  +  (miAXL-l)/  '/'t€iAX-l  /  J 


?noo  . 
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3.b  General  Planform  Geometry  Gescriptiori 


In  order  to  discuss  the  input  procedure  we  moat  first  define  a 
reference  system.  Let  us  define  a  cartesian  coordinate  system 
(Xl,  Yl,  Zl)  such  that  the  rotor  blade  axis  of  rotation  is  the  7.1 
axis.  The  projection  of  the  blade  span  axis  into  the  (XI,  Yl)  plane 
will  lie  alonR  the  Yl  ajcis  and  the  XI  axis  w:' 11  f^ene-rally  be  directed 
toward  the  pointed  trailing'  edge.  The  blade  must  have  a  non  zero 
hub  radius.  Now  with  the  rotor  blade  in  its  lifting  configuration 
within  this  reference  system  one  must  be  able  to  define  a  system  of 
surface  grid  point  coordinates  as  follows: 

A.  inscribe  a  series  of  non  intersecting  lines  (not  necessarily 
straight  lines)  on  the  rotor  blade  uppei  acrface  Joining  a 
point  on  the  root  section  to  a  point  on  the  tip  section. 

Two  of  these  lines  must  lie  along  the  leading  edge  and  traiJing 
edge.  Index  these  lines  consecutively  as  N  =  1 .'il-lAXU  where 
N  =  1  is  the  line  along  the  leading  edge  and  N  ^  NilAXU 
la  the  line  along  the  trailing  edge.  HllAXU  must  he  no  greater 
than  hi. 

.1.  Repeat  (A)  above  for  the  lower  surface  where  n  jw  N  =  i, 

NtvlAXL  hi. 

C.  Inscribe  a  aeries  of  non  intersecting  lines  on  the  rotor  blade 
upper  and  lower  surface  starting  at  a  point  on  the  leadin/: 
edge  going  to  the  trailing  edge  along  the  upper  surface  and 
returning  to  the  leading  edge  point  along  the  lower  surface. 
Identify  these  lines  as  M  =  1,  MMAX  -  k9  where  M  =  1  lies 
along  the  root  and  M  =  lies  along  the  tip. 

U.  The  intersection  of  a  line  of  constant  M  with  a  line  of 
constant  N  will  now  define  a  control  grid  point  indexed  as 
(N,M).  Four  lines  Joining  grid  points  (N,M),  (N  ,7+1 ),  f  .,+i  , 

(N+1,M)  will  define  a  control  surface  such  that  the  sum  of 
the  control  surface  areas  equals  the  wetted  area  of  the  blade. 
The  number  of  control  surface  areas  must  not  exceed  2000,  that 
is  /(NMAXU-1)  +  (NMAXI^l)  /  x  (M?iAX-l)^-  2000. 

L.  .<ow  for  each  grid  point  (N,M)  the  blade  sirface  coordinate:;  in 
terras  of  (Xl,Yl,Zl)  must  be  defined.  These  lengths  cire  to 
be  non  dimensional! zed  on  the  rotor  radius. 
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Having  determined  this  information,  the  data  input  is  card 
punched  and  arranged  as  follows: 


1.  Punch  on  one  card  the  value  of  !JMAXU,  NliAXL  and  R’^lAX  accordinp 
to  the  FORIiAT  (3110). 

P.  Puncli  a  set  of  carcfe  such  that  the  complete  set  of  upper 
surface  grid  coordinates  at  all  (N,M)  grid  points  is  given 
on  these  cords  as  H,M,X1U(N.M)  ,  YlU(ri,M),  ZlU(rJ,M) 
according  to  the  FORIIAT  (215 ,3E20. 6 ) . 

3.  Repeat  (2)  above  for  the  complete  set  of  lower  surface  grid 
coordinates . 


I*. a  Rotor  Free  Ctream  and  Wake  Analysis  Parameters 


'^he  following  data  is  to  be  inputted  in  the  NAfIKLII ■I'/DATAl/  data 
statement ; 


LAMB])A  =  A. a  (Default  value  is  O.C)) 

where  L/AMdUA  is  the  free  stream  rotor  inflow  ratio  and  is  given 
positively  by  the  rotor  climb  speed  divided  by  the  rotor  tip  s])eea. 


WAV  =  A. a  (Default  value  is  O.j) 

where  WAV  is  the  average  induced  downwash  at  ttie  trailing,  eig.e 
which  may  be  determined  using  actuator  disk  theory  assuming  a 
given  thrust.  This  parameter  along  with  LAMBDA  is  used  to 
prescribe  the  wake  geometry. 


JVAI’I  =  A. a  (Default  value  is  ?0.0°; 

where  D'Wfl  ir,  the  angle  increment  in  degrees  which  determines  the 
initial  wake  elemental  areas  used  for  calculating  near  wake 
influences . 


DWfP  =  A. a  (Default  value  is  D'/i'l ) 

wiiere  D'Wi’2  is  the  angle  increment  in  degrees  which  determines  the 
wake  elemental  areas  used  for  calculating  far  wake  influences. 


ANGLE  =»  A.a 


(Default  value  is  100.0°) 


where  AJJGIJi  is  the  angular  measure  in  degrees  from  trie  span 
axis  beyond  which  point  the  wake  calculations  are  to  be  done 
in  increments  of  DV/T?. 


'i'UHNS  =  A.a  (Default  value  i;,  1.0) 

where  TURNC  is  the  number  of  wake  revolutions  from  the  span  axis 
beyond  which  the  wake  calculations  are  to  cease. 


INFLU  =  K  (Default  value  is  0) 

where  INFLU  is  an  option  which  allows  the  program  to  cease  wake 
calculations  at  an  angle  comparable  to  TURNU  above  if  the  wak.c 
elemental  influence  coefficients  becomes  less  than  a  prescribed 
amount  described  in  AIHFLU  belo’  .  Select  K  as  follows: 

K  =  1  to  effect  this  option 
K  =  0  to  negate  this  option 

AINFLU  =  A.a  (Default  value  is  O.Oi) 

where  AiNFLU  is  a  parameter  defined  as  the  absolute  value  ratio 
of  the  total  influence  at  a  point  of  one  wake  element  to  the 
total  influence  of  another  wake  element.  The  point  at  which 
the  influences  are  compared  is  the  centroid  of  the  first  span 
station  upper  leading  edge  element.  The  program  computes  the 
influence  for  all  elements  of  wake  span  station  one  comparing 
this  always  to  the  influence  of  the  first  wake  element.  When 
this  ratio  becomes  less  than  AINFLU  the  parameter  TURNU  is 
redefined  to  effect  the  termination  of  calculations  at  this 
point  for  all  span  stations. 


lUAFi;  =  K  (Default  value  ir  U50000) 

where  IGAFE  is  a  parameter  which  indicates  the  total  number 
of  wake  elemental  influence  calculations  that  are  to  be  performed. 
Internally  the  program  computes  from  this  the  meocimum  number  of 
wake  elements  trailing  any  one  blade  at  one  span  station  that 
will  approximate  this  ISAFE  number  of  calculations  thus  it  ensures 
all  span  wake  segments  have  the  same  number  of  spiral  wake  elements. 


NLIFT  =  K 


(Default  value  is  O) 
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wh'  e  flLIPr  ia  a  paxameter  which  indicates  whether  the  program 
is  )  include  a  waJce  analysis  (Lift  Case)  or  is  not  to  include 
a  wake  analysis  (No  Lift  Case).  Select  K  as  follows: 

K  =  1  NO  LIFT  OPTION 

K  =  0  LIPr  OPriON 

U4AX  =  K  (Default  value  is  ?) 

where  LIIAX  is  the  number  of  equally  spaced,  identical  rotor 
blades. 

Note  that  j  n  the  above  data  there  are  three  input  parameters 
(turns,  INFLU.  IGAFE)  which  may  be  used  to  terminate  the  wake  analysis 
calculations.  They  are  in  fact  independent  of  each  other, 

l4,b  Planar  Free  Stream  and  Woke  Analysis  Parameter 

The  following  data  is  to  be  inputted  in  the  NAI'tELIST/UAVAft/  data 
statement : 

WAV  =  A.a  (Default  value  is  0.0) 

where  WAV  is  the  average  induced  downwash  at  the  trailing  edge 
nondimensionalized  by  the  free  stream  velocity.  This  velocity 
is  used  to  wash  the  wake  below  the  plane  of  the  free  stream 
velocity. 

NLIFT  =  K 

where  NLIFT  is  as  defined  in  the  preceeding  section. 

Please  note  that  the  planar  wing  analysis  option  of  this  program 
has  not  at  the  present  time  been  thoroughly  checked.  It  is  included 
here  only  to  document  the  program.  The  planar  wing  analysis  option 


is  obtained  by  altering  the  rotor  blade  analysis  routine  and  we 


l6') 


acknowledge  the  fact  that  this  approach 


is  very  inefficient. 


D.P.I4  Prograjn  Output 


The  main  output  from  this  program  is  as  follows: 


A.  Printed  output  which  essentially  prints  out  the  data  that 
was  read  into  the  program. 

B.  Printed  output  which  essentially  is  a  program  monitor. 

This  output  indicates  the  cuinmulative  time  at  which  the 
major  subroutines  were  entered  and  whether  they  terminated 
normally.  Additional  output  indicating  specific  itpm;i  are 
included  in  the  CBIjAPe:  and  CWAKE  related  monitor  output 
statements.  This  output  is  self  explanatory. 

C.  Data  output  stored  on  external  files  ITl  and  IT2.  Unlike 
items  (a)  and  (B)  above  this  output  is  required  for  the 
continuation  of  the  overall  program.  We  will  describe  the 
detailed  data  on  these  two  files  when  we  discuss  the 
"MAIN-PniNT"  program  whose  function  is  to  retrieve  specific 
data  from  these  files  and  print  on  paper  this  data. 


In  addition  to  the  above  listed  output  additional  output  rn;iy  be 
selected.  We  have  already  discussed  the  nature  of  this  output  in  tne 
section  describing  the  input  procedures  and  program  options  available, 
he  ujtput  as  printed  is  self-explanatory. 
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D.3  MAIN  -  CREATE  SOLUTION  FILE  PROGRAM 

In  this  section  we  will  discuss  the  program  MAIN  -  CREATE 
SOliUTION  FILE.  This  program  is  rather  straightforward  with  it;; 
function  being  to  operate  on  the  coefficient  data  written  on  external 
files  ITl  and  ITS  by  the  MAIN  -  INFLUENCE  CALCULATION  program,  in 
order  to  generate  the  set  of  normal  influence  coefficients  and 
normal  free  stream  velocities.  This  set  of  coefficients  and  velo¬ 
cities  are  then  stored  on  external  file  IT3  or  IT3  and  ITh  in  a  form 
compatable  for  the  solution  programs  which  will  be  discussed  in  the 
next  section. 

D.  .1  Guramary  Discussion  of  the  Subroutines 

•ubroutlnes  SETUP2  (ITl,  IT2,  IT3) 

The  function  of  this  subroutine  is  to  generate  from  the  data 
given  on  externeLl  files  ITl  and  IT2  the  set  of  normal  influence 
coefficients  expressed  as  an  IMAX  square  matrix  [AD  and  to  generate 
the  set  of  IMAX  normal  free  stream  velocities  expressed  as  the  negative 
of  the  IMAX  matrix  IBIsuch  that  the  linear  set  of  IMAX  algebraic 
equations  to  be  solved  is  of  the  form 

[A1  |x|  ^  IBI 

This  data  is  then  written  on  external  file  IT3  in  double  precision  as 
follows  where  each  line  represents  a  separate  write  statement: 


ITi 

IMAX 

B(l),  B(2),  (B3) 

A(l,l),  A(1.2)  . 


B(IMAX) 
A(l,  IMAX) 


AdMAX.l) . A(IMAX,II-IAX) 

This  data  file  is  compatable  with  the  exatt  solution  program  discussed 

later . 


Subroutine  SETUPS  (IT1,IT2,IT3,IT**) 

The  fiinction  of  this  subroutine  is  the  same  as  SETUP2  to  the  point 
where  the  algebraic  equations  are  defined  by 

[AJ  JXI  =  IB j 


This  subroutine  further  operates  on  the  coefficient  such  that  the  set 
of  equations  may  be  written  as 


l-J  J 


which  in  matrix  form  becomes 


/A  J 


{A  J  /Xj  f  iOj 


where 


O- 


^7 


4  -* 


/?yj 


// 


// 


// 


17P 


This  subroutine  then  writes  on  external  file  IT3  the  [Cj  matrix 
as  shown  below  where  each  line  represents  a  separate  write  statement: 


4/ 

4/ 

- 

.  • 

. 

. 
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On  external  file  ITU  is  written  the  colvunn  matrix  |dI  . 

This  data  as  written  on  files  IT3  and  ITU  is  compatable  with  the 
iteration  program  discussed  later. 

Subroutine  VWl 

See  dicussion  in  section  D.2.1. 

D.3.2  Program  Operation 

The  program  operation  requires  three  or  foiur  external  files 
depending  on  whether  the  solution  will  be  obtained  using  the  exact 
elimination  or  the  iterative  solution  method  respectively.  If  three 
files  are  required  then  thv^  foiurth  file  must  be  defined  as  a  dummy 
file.  Since  the  main  program  simply  calls  on  either  one  of  two  sub- 
routines  no  flow  diagrai^  will  be  presented. 

D . 3 . 3  Program  Input 

The  input  required  consists  of  the  following  parameters  given 
as  part  of  a  NAMELIST/SOLUF/  data  statement: 

METH  =  K 

where  METH  is  a  control  parameter  which  indicates  the  proper 
solution  file  tc  be  created.  Select  K  as  follows: 

K  =  1  Elimination  Solution  File 
K  =  2  Iterative  Solution  file 


|;iijyji^,,J'tti.WV-  U’ 
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ITl  =  K  (Default  value  is  l) 

whore  K  la  the  reference  number  assigned  to  the  comparable 
external  file  ITl  of  the  Main  -  Influence  Coefficients  program. 
No  changes  will  be  made  to  the  file  contents. 


IT2  =  K  (Default  value  is  2) 

where  K  is  the  reference  number  assigned  to  the  comparable  external 
file  IT2  of  the  Main  -  Influence  Coefficient  program,  ilo  changes 
will  be  made  to  the  file  contents. 


IT3  =  K  (Default  value  is  3) 

where  K  is  the  reference  niunber  assigned  to  the  input  file 
created. 


ITl*  5s  K  (Default  value  is  U) 

where  K  is  the  reference  number  assigned  to  the  additional  solution 
file  required  if  METH  =  0  option  is  used.  If  Ml^TH  =  1  option  is 
chosen  input  the  same  reference  niunber  for  ITl*  as  inputted  for  IT3. 

D.3.1*  Progreun  Output 

The  output  of  this  program  consists  of  the  following: 

1)  A  printed  statement  indicating  whether  the  elimination  or 
iterative  solution  file  was  created. 

2)  Data  written  on  external  file  IT3  if  the  elimination  solution 
file  was  created  or  data  written  on  external  file  IT3 

and  ITl*  if  the  iterative  solution  file  was  created. 
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D,U  MAIU  -  SOLVE  BY  ELIMINATION  PROGRAM 

In  this  section  we  will  discuss  the  program  tiAIN  -  SOLVE  BY 
ELIMINATION.  The  function  of  this  program  is  to  solve  a  simultaneous 
set  of  equations  using  the  Gaussian  elimination  method  with  pivotal 
condensation  as  described  in  Section  3.^.  This  program  is  not 
competitive  with  routines  not  requiring  peripheral  storage  and  it 
is  suggestv'fd  that  this  method  be  used  only  for  simultaneous  sets  of 
equations  having  less  than  350  unknowns  because  of  time  considerations. 

D.U.l  Summary  Discussion  of  the  Subroutines 
Subroutine  SSIMQ 

The  function  of  this  subroutine  is  to  accept  the  coefficients  as 
written  on  an  external  file  in  the  manner  described  in  the  discussion 
of  subroutine  SErUP2  (Section  D.3.1)  and  operate  on  this  matrix  of 
coefficients  in  order  to  reduce  it  to  an  effective  diagonal  matrix 
whose  lower  half  is  composed  of  zero  elements.  This  is  the  pivot  eli¬ 
mination  part  of  the  solution  routine. 

Subroutine  BSUB 

The  function  of  this  suoroutine  is  to  operate  further  on  the 
diagonal  matrix  derived  from  SSIMQ  and  by  the  back  substitution  method 
obtain  the  unknowns  to  the  simultaneous  set  of  equations. 

Subroutine  HEAD 

Refer  to  Section  D.2.1. 
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D.U.2  Program  Operation 

The  operation  of  this  program  is  relatively  straightforward. 

It  does  require  four  or  three  external  data  files  depending  on  whether 
the  Initial  coefficients  are  to  be  saved  ’jntll  the  solution  is  completed 
or  not  saved  respectively.  The  InitlsLl  data  must  be  giv-ii  on  file 
ITF  and  the  final  results  are  written  on  file  ITF.  Files  labeled 
ITA,  and  ITC  should  be  high  speed  data  files  as  they  are  actively  used 
in  the  elimination  phase  of  the  program.  If  the  initial  data  is  n-^t 
to  be  kept,  file  ITF  and  file  ITA  may  be  the  same  file.  File  ITB  is 
used  to  record  the  diagoneil  matrix  smd  subsequently  is  used  along  with 
file  ITF  in  the  back  substitution  phase  of  the  problem.  The  speed  of 
transferring  data  from  files  ITB  and  ITF  into  and  out  of  main  core 
storage  is  not  as  critical  as  for  ITA  and  ITC. 

The  program  does  contain  a  series  of  options  which  allow  the 
solution  to  proceed  in  a  series  of  executable  runs  or  in  one  single 
executable  run.  These  options  are  described  in  the  next  section  but 
briefly  we  may  list  the  options  here  as  follows: 

A.  The  program  may  proceed  from  data  entry  through  the  elimination 
phase,  through  the  back  substitution  phase  to  the  recording  of 
the  solution  value. 

B.  The  program  may  proceed  from  data  entry  into  the  elimination 
phase  and  terminate  after  a  given  amount  of  time  has  elapsed, 
or  after  the  elimination  phase  ■'s  completed. 

C.  The  program  may  be  restarted  in  the  elimination  phase  and  terminated 
as  in  (E). 

D.  The  program  may  initiate  execution  in  the  back  substitution  phase. 
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The  solution  values  as  written  on  file  ITF  are  in  a  form  compatable 
for  input  into  the  program  MAPI  -  MOVE  oGLUTION  FILE  which  transfers 
the  answers  to  one  of  the  two  permanent  files . 


n.4.3  Program  Input 

The  following  parameters  are  to  be  inputted  as  part  of  the 
NAMELILT/GELIM/  data  statement: 


IP  =  K 


where  IP  is  aui  option  control  parameter  such  that  for: 

K  =  0  The  prograi.i  wi..l  solve  the  system  of  equations  in  one 
execution  step. 

h  =  1  The  program  will  begin  execution  in  the  back  substitution 
routine  assumlnr  the  elimination  phase  has  been  completed 
and  the  necessary  data  is  on  file  ITB  described  later. 

h  The  program  will  begin  execution  as  for  K.  =  0  above 

but  will  terminate  after  the  elimination  phase  is 
completed. 

K  a  3  The  program  will  restart  In  the  elimination  routine  at 
elimination  row  (M)  and  pivot  row  (N)  described  later 
in  the  TMAX  discussion. 


TMAX  =  A. a 


(Default  value  is  1200.0) 


where  TMAX  is  a  parameter  which  is  an  input  estimate  of  the  time 
in  seconds  required  for  the  elimination  phase  of  the  program. 

If  the  elimination  phsLse  has  not  been  totally  completed  by  this 
time  the  program  will  terminate  the  elimination  phase  after  setting 
up  files  ITA  and  ITB  in  a  manner  compatable  for  the  restart  option 
selected  by  IP  =  3  above.  All  restart  inforiiation  including  the 
resteurt  elimination  row  (M)  and  pivot  row  (N)  are  written  for 
files  ITA  and  ITB. 
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i;<AX  =  K  (Default  value  is  0) 

where  IMAX  is  a  parameter  Indicating  the  number  of  unknowns  in 
the  system  of  equations. 


DUMMY  =  K 

where  K  must  be  inputted  as  69^6  if  TMAX  above  is  greater  than  or 
equal  to  1200. 


ITA  =  K. 


(Default  value  is  K  =  1) 


where  K  is  the  reference  number  assigned  to  external  file  ITA 
which  is  a  data  file  active  in  the  eliminatior  routine  and 
contains  data  necessary  for  the  re8teu*t  option. 


ITB  =  K  (Default  value  is  K  =  ?) 

where  K  is  the  reference  number  assigned  to  external  file  ITb 
which  is  a  data  file  active  in  the  elimination  routine  and  contains 
data  necessary  for  the  restart  option  and  the  back  substitution 
routine. 


ITC  =  K  (Default  value  is  K  =  3) 

where  K  is  the  reference  number  assigned  to  external  file  I'l'C  which 
is  a  temporary  data  file  active  in  the  elimination  phase.  The 
data  written  on  it  is  not  necessary  for  the  restart  option  or  back 
substitution  routine. 


ITF  =  K  (Default  value  is  ITA) 

where  K  is  the  reference  number  assigned  to  external  file  ITF  which 
is  the  initial  data  input  file  created  by  program  "MAIN-CRKATi. 
SOLUTION  FILE"..  This  is  adso  the  file  on  which  the  final  solution 
values  are  written  replacing  the  column  input  data  record.  File 
ITA  may  be  taken  to  be  the  same  as  file  ITF  but  in  this  case  the 
data  is  destroyed  before  the  program  has  obtained  the  solution 
values . 
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D .  .  1*  Program  Output 

The  main  output  of  this  program  conalata  of  the  following: 

A.  A  printed  statement  indicating  whether  a  successful  solution 
or  singular  solution  resulted  or  a  statement  indicating  where 
the  solution  was  terminated. 

U.  A  printed  listing  of  the  solution  vedues  which  are  also 

written  on  File  ITF  in  a  form  compatable  for  data-input  into 
the  program  MAIN  -  TRAT'cffiR  SOLUTION  FILE. 

If  the  elimination  phase  is  to  be  rim  in  stages ,  files  ITA  and 
ITB  should  be  considered  as  output  both  of  which  are  needed  as  the  data 
set  input  for  the  restart  option. 

If  the  program  terminates  after  the  elimination  phase  then  file 
ITB  should  be  considered  as  output  which  becomes  the  data  set  input  for 
the  bach  substitution  phase. 
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D.5  MAIN  -  GOLVE  BY  ITERATION  PROGRAM 

In  this  section  we  will  discuss  the  program  MAIN  -  SOLVE  BY 
ITERATION.  The  function  of  this  program  Is  to  solve  a  simultaneous 
set  of  equations  using  the  simultaneous  displacement  Iterative  method 
described  In  section  3.6. 

D.5.1  Summary  Discussion  of  the  Subroutines 

Subroutine  ITER ( FB ,FC , IMAX  ,K ,KMAX  .DELTA  ,TMAX  , JR  ) 

The  function  of  this  subroutine  Is  to  perform  the  actual  Iterative 
procedure.  The  IMAX  set  of  coefficients  are  given  on  external  files 
labeled  FB  and  FC  which  are  described  in  section  D.3.1  in  the  dis¬ 
cussion  of  SETUP3.  The  parameter  K  indicates  the  iteration  number. 

Delta  is  a  test  parameter  which  dictates  the  greatest  difference  that 
all  present  solution  values  may  differ  from  their  previous  Iterative 
values  in  order  for  the  Iteration  scheme  to  be  completed.  KMAX  and 
TMAX  are  parameters  which  may  terminate  the  iterative  scheme  after 
KMAX  iterations  or  after  TMAX  elapsed  seconds.  The  parameter  IR  is  a 
return  code  indicating  on  what  basis  the  subroutine  execution  was 
tex*mlnated. 

Subroutine  DIFF  (FB.FC ,IMAX ,K) 

The  function  of  this  subroutine  is  to  compute  the  value  of  the 
solution  vector  using  the  coefficients  and  compare  this  value  to  the 
iterative  solution  vector  and  print  out  the  results. 


Subroutine  HEAD(K) 

Refer  to  section  D.2.1. 

D . 5 • 2  Program  Operat 1 on 

The  program  operation  requires  two  external  data  files.  These 
files  are  the  solution  files  created  in  the  MAIN  -  CREIATE  SOLUTION 
FILE  program  emd  are  discussed  more  specifically  in  the  SETUP3 
discussion  of  section  D.3-1* 

The  program  contedns  a  number  of  options  which  permit  the  program 
to  terminat<"  when  a  certain  convergence  criteria  is  established,  or 
when  a  given  number  of  iterations  have  been  performed  or  when  a  given 
amount  of  time  has  elapsed.  The  program  prior  to  termination  writes 
out  on  one  of  the  external  data  files  the  latest  iterative  solution 
values  thus  allowing  the  program  to  be  executed  again  using  these 
last  values  as  the  initial  iterative  values  for  another  series  of 
Iterations.  The  solution  values  are  written  on  file  FC  immediately 
after  the  data  records  written  on  file  FC  by  GETUP3  and  they  are 
written  in  a  form  compatable  for  input  into  the  progroi.i  .’lAI.!  -  TRAiiGFLR 
GOLUTIO;!  FILE  which  transfers  the  solutions  to  one  of  the  two  permanent 
files. 

U.5.3  Program  Input 

The  following  parameters  are  to  be  inputted  as  part  of  the 
NAMELIGT/SITER/  data  statement: 
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FB  =  N  (Default  value  is  N  =  1) 

where  N  is  the  reference  number  indicating  the  external  file 
on  which  the  square  coefficient  matrix  is  written.  This 
data  file  is  data  file  IT3  created  in  the  MAIN  -  CREATE 
SOLUTION  FILE  program  if  option  METH«0  was  executed.  (Refer 
to  section  D.3) 


FC  *  N  (Default  value  is  N  =  2) 

where  K  is  the  reference  number  indicating  the  external  file 
on  which  the  column  matrix  is  written.  This  data  file  is 
data  file  IT**  created  in  the  MAIN  -  CREATE  SOLUTION  FILE 
program  if  option  METH=0  was  executed.  (Refer  to  section  D.3). 
Ihiring  execution  the  present  iteration  solution  values  are 
stored  on  this  file  immediately  after  the  column  matrix  record 
destroying  any  other  previously  written  Iterative  values. 


IMAX  =  N  (Default  value  is  N  =  0) 

where  N  is  the  number  of  unknowns. 


K  =  N  (Default  value  is  N  =  0) 

where  the  parameter  K  dictates  whether  previous  iterative 
solution  values  are  stored  on  file  FC  and  are  to  be  used  as 
starting  iterative  values.  Select  N  as  follows 

N  =  1  Indicates  previous  solution  values  are  stored 
on  file  FC  and  are  to  be  used  as  starting 
iterative  solution  values. 

N  =  0  Indicates  no  previous  solution  values  are  stored 

on  file  FC  and  the  program  will  assume  the  starting 
iterative  solution  values  to  be  all  zero. 


iCMAX  =  N  (Default  value  is  N=K+10) 

where  the  value  of  W  indicates  the  number  of  iterations  to  be 
performed  subject  to  DELTA  and  TMAX  constraints  described 
below. 
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TMAX  =  A. a  (Default  value  is  120.0) 

where  A. a  Is  the  cummulative  time  in  seconds  after  which  the 
program  is  to  terminate  the  iteration  scheme. 


DELTA  »  A. a  (Default  value  la  O.Ol) 

where  DELTA  is  a  convergence  criteria  such  that  if  all 
calculated  values  of  the  solution  vector  less  the  Iterative 
solution  vector  values  divided  by  the  calculated  solution 
vector  values  are  less  than  DELTA  in  the  absolute  sense, 
then  the  iteration  scheme  is  terminated. 


D.5.^  Program  Output 

The  printed  output  of  the  program  consists  of  the  following: 

1.  A  printed  line  of  output  after  each  iteration  which  indicates 
the  Iteration  number,  the  number  of  present  solution  values 
whose  change  from  the  previous  solution  values  are  less 

them  the  prescribed  convergence  parameter,  the  cummulative 
time  and  the  actual  value  of  the  first  number  in  the  solution 
vector  and  its  change  from  its  previous  value. 

2.  A  printed  line  indicating  on  what  basis  the  program 
terminated. 

3.  After  the  last  iteration  the  entire  solution  vector  and 
the  difference  vector  (present  solution  less  previous 
solution)  cure  printed. 

In  addition  to  the  abive  the  solution  vector  is  also  written 
on  file  FC  and  this  output  should  be  considered  as  input  for  successive 
iteration  program  executions  or  for  input  to  the  MAIN  -  TRANSFER 
SOLUTION  FILE  program. 
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D.6  MAIN  -  TRANSFER  SOLUTION  FILE  PROGK\M 

The  function  of  the  MAIN  -  TRANSFJB  SOLUTION  FILE  propram  Is  to 
read  the  solution  values  written  on  an  external  file  by  either  the 
elimination  or  iteration  scheme  program  and  to  write  these  same 
answers  on  one  of  the  two  permanent  files  so  that  the  two  permsuient 
files  now  have  recorded  on  them  eDl  the  information  necessary  to 
compute  the  various  velocities,  pressures  amd  forces. 

The  program  requires  two  files  for  its  operation  and  is  rather 
straightforward.  The  input  required  for  this  program  includes  the 
specification  of  the  following  two  parameters  as  part  of  the 
NAMELIST /MOVEC/  data  statement; 


FA  s«  K  (Default  value  is  K  =  1) 

where  K  is  the  reference  number  assigned  to  the  permanent  file 
FA  which  was  initially  created  in  the  MAIN  -  INFLUENCE  CALCU¬ 
LATIONS  program  and  was  there  designated  file  ITl. 


FC  =  K  (Default  value  is  K  =  E) 

where  K  is  the  reference  number  assigned  to  the  file  FC  which  has 
written  on  it  the  solution  vector.  If  the  solution  was  obtained 
using  the  elimination  method  then  file  FC  is  identical  to  tne 
elimination  file  ITF.  If  the  solution  was  obtained  using  the 
iterative  scheme  then  file  FC  is  identical  to  the  iterative 
file  FC. 


The  printed  output  of  this  progreun  is  simply  a  statement  indi¬ 
cating  successful  completion  of  the  program. 
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D.7  MAIN  -  VELOCITY  CALCULATIONS  PROGRAM 


In  this  section  we  will  discuss  the  MAIN  -  VELOCITY  CALCULATIONS 
prot’ireun.  Tlie  function  of  this  prof^ram  is  to  compute  on  each  blade 
element  1=1,  IMAX,  the  various  velocities,  pressure  coefficients 
and  forces  and  to  print  these  results  as  well  as  a  summary  of  the 
geometry  for  each  of  these  blade  elements.  The  program  also  computes 
and  prints  the  rotor  thrust,  torque  and  in  plane  force  coefficients. 
The  subroutines  called  are  WTl  and  HEAD.  These  have  previously  been 
discussed  in  section  D.2.1  and  will  not  be  repeated  here. 

The  program  requires  for  its  operation  the  two  external  permanent 
files  initially  created  in  the  MAIN  -  INFLUENCE  CALCULATIONS  program 
emd  a  temporary  file  to  be  used  during  the  execution  of  this  program. 

D . 7 . 1  Program  Input 

The  program  receives  all  of  its  data^ input  from  the  permanent 
file.  The  card  data  input  Includes  specification  of  the  following 
parameters  as  paurt  of  the  NAMiLLIST/DATAV/  data  statement: 

ITl  =  K  (Default  value  is  K  =  1) 

where  K  is  the  reference  number  assigned  to  the  permanent 
external  file  ITl.  This  file  is  to  be  identical  to  file  ITl 
created  in  the  MAIN  -  INFLUENCE  CALCULATIONC  program. 


IT2  =  K  (Default  value  is  K  =  2) 

where  K  is  the  reference  number  assigned  to  the  pemnanent  external 
file  IT2.  This  file  is  to  be  identical  to  file  IT2  created 
in  the  MAIN  -  INFLUENCE  CALCULATIONS  program. 
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IT3  =  K 

where  K  is  the  reference  nvmber  assigned  to  the  temporary 
external  file  IT3.  This  file  is  used  during  execution  of  the 
program  and  upon  completion  contains  no  useful  information. 


D.T.?  Program  Output 

The  program  prints  out  the  following  summary  results  using  the 
indicated  symbols  for  each  blade  surface  element  1=1,  IMAX; 


1.  Element  index  location  in  terms  of  I  and  in  termls  of  indices 
(rJ,M)  where  N  is  the  chord  station  and  M  is  the  span  station, 

2.  Centroid  location  in  terms  of  the  blade  fixed  coordinate  system. 

The  symbols  used  is  RCG  where  RCG  =  (XCG^  +  YCG^)^, 

3.  The  length  of  the  curve  connected  by  the  centroids,  starting  from 
the  leading  edge,  following  the  chordwise  direction  and  ending 

at  the  trailing  edge.  The  symbol  used  is  CHODTL. 

U.  Nine  components  of  the  transformation  matrix  between  the  elemental 
surface  coordinate  system  and  the  blade  fixed  coordinate  system. 
The  symbols  used  are  TPIX,  TPIY,  TPIZ,  TP2X,  TP2y,  W,  NY,  NZ. 

5.  The  chordwise  and  spanwlse  surface  velocities  and  their  resultant 
at  the  centroid  of  the  element  with  respect  to  tlie  blade  fixed 
coordinate  system.  The  symbols  used  are  VTl',  VT2',  and  VT'. 

6.  The  free  stream  velocity  and  its  ccmponents  at  che  centroid  of 
the  element  with  respect  to  the  blade  fixed  coordinate  system. 

The  symbols  used  are  respectively  VINF',  VINFX',  VINFY ' ,  and  VINFZ'. 

7.  The  components  of  the  free  stream  velocity  in  chordwise  and 
spanwise  directions  at  the  centroid  of  the  element  with  respect 
to  the  blade  fixed  coordinate  system.  The  symbols  used  are 
VINFTl'  and  VINFT2'. 

8.  The  derivatives  of  the  velocity  potential  in  both  chordwise  and 
spanwise  directions.  The  symbols  used  are  DPHIDC  and  DPHTDY. 

9.  The  pressure  coefficient,  CPI,  nondimensionalized  on  the  basis 
of  the  tip  speed  velocity. 

10.  The  pressure  coefficient,  CP2 ,  nondimensionalized  on  the  basis 
of  the  local  free  stream  velocity. 


11.  The  ratio  of  the  surface  velocity  to  the  free  stream  velocity 
(VT/VINF). 

12.  The  double-*,  strengths,  MU. 

13.  The  velocity  potentials,  PHI. 

1*4.  The  angle  that  the  surface  velocity  raetkes  with  the  chordwise 
direction,  THETA. 

15.  The  elemental  force  components  (FXE,  FYF,  and  F7.F)  in  terms  of 
the  blade  fixed  coordinate  system. 

16.  The  elemental  torque  (QZE)  about  the  axis  of  rotation. 

All  lengths  and  velocities  above  are  nondimensionalized  on  the  basis 
of  the  rotor  radius  and  tip  spee <  respectively. 

In  addition  to  the  above  a  summary  of  the  rotor  aerodynamic 

coefficients  are  printed.  This  output  is  self  explanatory.  Also 

the  differences  between  the  pressure  coefficients  of  upper  and 

lower  surfaces  are  calculated  and  printed. 


U.8  MAIN  -  PPINT  PROGRAM 


In  this  section  we  will  discuss  the  MAIN-  PRINT  program.  The 
function  of  this  program  Is  to  selectively  read  data  stored  on  the 
two  permaiient  data  set  files  auid  print  this  data.  This  is  aji  auxiliary 
program  and  is  not  used  in  the  normal  execution  of  the  overall  program. 
It  was  designed  initially  to  be  used  in  checking  out  the  program 
operation.  We  will  describe  here  how  this  program  may  be  used  for  the 
purpose  of  printing  out  the  net  blade  influence  or  wake  influence 
coefficients  on  some  blade  control  element  I,  and  for  printing  out  a 
complete  description  of  the  elemental  control  surfaces.  The  program 
does  require  for  its  operation  the  two  permanent  external  flics 
created  in  the  MAIN  -  INFLUENCE  CALCULATIONS  program. 

D.8.1  Program  Input 

The  card  input  data  consists  of  the  following  parameters  inputted 
as  peu't  of  the  NAMELIST/DATAP/  data  statement: 

ITl  =  K  (Default  value  is  K  =  1) 

where  K  is  the  reference  number  assigned  to  file  ITl  which  is 
identical  to  file  ITl  created  during  the  MAIN  -  INFLUENCE  CALCU¬ 
LATIONS  program. 


IT2  =  K  (Default  value  is  K  =  2) 

where  K  is  the  reference  number  assigned  to  file  IT2  which  is 
identical  to  file  IT2  created  during  the  MAIN  -  INFLUENCE  CALCU¬ 
LATIONS  progreim. 
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ISURF  «  K  (Default  value  is  K  =  0) 

where  K  is  an  integer  indicating  whether  the  control  3urfac^: 
geoDietry  data  is  to  be  printed.  Select  K  as  follows: 

K  »  1  Surface  data  printed 
K  ■  0  '  Surface  data  not  printed 


Ml  *  Kj^,  K2,  .  .  •  Kp.  .  .  IC20  (Default  value  is  Kp  =  0) 

where  Kp  >  Kp_i  and  Kp  is  the  integer  spanwlse  wake  station  for 
which  the  wake  influence  coefficients  on  all  blade  elements 
1=1,  IMAX  are  to  be  printed.  Select  K  as  follows: 

=  0  No  wake  Influence  coefficient  piinting. 

^2^  =  9999  All  spanwise  wake  influence  coefficients  are  printed. 

=  K  The  spanwise  waike  influence  coefficients  are 
printed. 


II  «  K]^,  K2,  .  .  .  Kp,  .  .  .  K20  (Default  value  is  Kp  =  O) 

where  Kp  >  Kp.^  and  Kn  is  the  integer  blade  control  surface  on 
which  all  blade  elemental  control  surfaces  Influence  coefficients 
are  to  be  printed.  Select  Kp  as  follows: 

Kj^  =  0  No  blade  Influence  coefficients  are  printed. 

=  9999  All  blade  influence  coefficients  are  printed. 

Kp  e  K  All  elemental  blade  influence  coefficients  on  blade 
control  surface  K  are  printed. 


D.8.2  Program  Output 

The  program  output  may  consist  of  the  following  depending  on  the 


options  selected: 
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A.  CONTROL  SURFACE  GEOMETRY  DATA 


This  conBlstB  of  the  geometry  data  which  descrlbeB  each  elemental 
control  surface.  It  Includes  the  following: 

1.  The  transformation  matrix  CaJ  relating  the  blade  fixed 
reference  system  to  the  elemental  coordinate  system.  The 
symbols  used  are 
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TPIX  TPIY  TPIZ 
TP2X  TP2Y  TP2Z 
TP3X  TP3Y  TP3Z 


2.  The  location  of  the  origin  of  the  elemental  coordinate  system 
with  respect  to  the  blade  fixed  coordinate  system  (XO,  YO,  ZO). 

3.  The  centroid  location  of  the  elemental  surface  with  respect 
to  both  the  blade  fixed  (XCG,  YCG,  ZCG)  and  elemental 
(XTCG,  YTCG)  reference  systems. 

>4.  The  area  of  the  element  (ACS). 

5.  The  corner  point  coordinates  of  the  elemental  surface  wit!, 
respect  to  the  elemental  reference  system,  (XTP,  YTP)^_j^ 

6.  The  ne.vative  of  the  free  stream  velocity  components  (VX, 

VY,  VZ)  at  the  elemental  centroid  location  with  respect  to 
the  blade  fixed  reference  system.  The  negative  of  the  normal 
free  stream  velocity  (VN)  is  also  given. 

B.  WAKE  INFLUENCE  COEFFICIENTS 

For  the  selected  spanwise  wake  station  the  influence  coefficients 
(CXW,  CYW,  CZW)  on  all  blade  control  surfaces  1=1,  IMAX  are  printed. 

These  influence  coefficients  are  referenced  to  the  blade  fixed  coordinate 
system. 


C.  BLADE  INFLUENCE  COEFFICIENTS 

For  the  selected  control  surface  the  influence  coefficients 
(CXB,  CYB,  CZB)  of  all  bleuie  elements  on  this  control  surface  are  printed. 
These  influence  coefficients  ore  referenced  to  the  blade  fixed  coordinate 
system. 


In  addition  to  the  above  output  those  parameters  of  common  statements 
LINDA  and  BRENDA  which  are  essentially  internal  control  parameters  are 
also  printed  out.  This  information  was  written  on  file  ITl  by  subroutine 
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Wri  in  the  MAIH  -  INFLUENCE  CALCULATIONS  program.  This  information 
is  shared  by  the  various  MAIN  programs  and  serves  as  the  continuity 
link  among  them. 


D.9  CONCLUDING  REMARKS 

We  have  attempted  in  the  previous  sections  to  describe  the  series 
of  computer  programs  and  their  operation  in  a  manner  more  detailed 
than  a  typical  "user's  manuaJ.".  We  have  chosen  to  describe  it  this 
way  in  order  to  document  not  only  the  program  but  the  philosophy  under 
which  the  program  systematically  performs  the  overall  problem.  With 
this  as  the  basic  program  we  hope  to  further  modify  the  various 
routines  in  order  to  optimize  the  program  in  terms  of  reducing  the 


computer  time  required. 


